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INTRODUCTION 


1.1  Summary 

Primary  activities  during  the  first  year  of  the  project  consisted  of: 

(i)  developing  a composites  processing  and  testing  laboratory  and  a neutron 
radiography  facility,  (ii)  conducting  research  in  accordance  with  the 
Statement  of  Work  given  in  Section  1.2,  (iii)  preparing  four  technical 
papers,  (iv)  various  interactions  of  the  faculty  with  the  technical  com- 
munity through  presentation  of  papers,  chairing  of  conference  sessions,  con- 
sulting, etc. 

Sections  2-4  describe  the  research  activities.  The  professional  per- 
sonnel associated  with  the  project  and  the  individual  activities  at  confer- 
ences are  given  in  Section  5.  Papers  completed  during  the  year  are  repro- 
duced in  the  Appendix. 

1.2  Statement  of  Work 

"a.  Investigate  the  effects  of  composite  cure  cycle  parameters  on  the 
mechanical  characteristics  of  structural  specimens  and  components  (see 
Section  2.): 

(1)  Investigate  the  development  of  new  innovative  in-process  monitor- 
ing devices  and  study  improvements  in  the  "Duomorph"  sensor. 

(2)  Conduct  experimental  measurements  on  both  neat  resin,  unidirec- 
tional laminates  and  cross-ply  laminates  to  establish  the  relationships 
at  various  elevated  temperature  and  relative  humidities  between  processing 
parameters  and  mechanical  strain  rate  behavior,  linear  and  identified 
nonlinear  viscoelastic  behavior,  dimensional  stability, and  fracture 
characterization. 

(3)  Investigate  the  influence  of  processing  parameters  on  Tg  at 
various  relative  humidities. 

b.  Investigate  the  influence  of  elevated  temperature  and  high  relative 
humidity  environments  on  composite  laminates  and  adhesively  bonded 
structural  components  (see  Section  3.): 

(1)  Investigate,  on  a theoretical  level,  the  formulations  of  a 
coupled  diffusion  theory  to  accurately  predict  the  environmental  degrad- 
ation of  composite  laminates  and  adhesive  mechanical  performance  which 
considers  the  influence  of  cure  cycle  parameters  being  Investigated  in 
a above. 


(2)  Predict  moisture  concentration  profiles  and  internal  stress/ 
strain  fields  within  a laminate  and  bonded  structure  as  well  as  in 
the  vicinity  of  exposed  edges.  These  profiles  will  consider  generic 
environmental  conditions  of  absorption  and  desorption. 

(3)  Compare  theoretical  predictions  with  available  experimental 
data  or  define  experiments  which  would  yield  necessary  data  for  this 
comparison. 

c.  Investigate  suitable  experimental  techniques  which  can  be  used  to 
nondestructively  examine  composite  laminate  initial  defects,  damage 
growth  and  environmental  degradation  characteristics  which  affect 
mechanical  performance.  Initially,  emphasis  will  be  given  to  neutron 
radiography  as  a suitable  technique  (see  Section  4.): 

(1)  Conduct  necessary  experiments  to  characterize  neutron  beam 
parameters  and  associated  parameters'  sensitivities  applicable  for 
nondestructively  examining  major  composite  laminate  mechanical  charac- 
teristics which  affect  strength,  dynamic  response  and  durability. 

(2)  Investigate  the  development  of  experimental  neutron  radiography 
techniques  to  measure  moisture  concentrations  in  composite  laminates  as 
well  as  spatial  concentration  profiles. 

(3)  Investigate  the  development  of  automated  procedures  for  the  anal- 
ysis of  data  from  neutron  radiography  experiments  which  accurately  and 
efficiently  characterize  and  identify  laminate  initial  defects,  mechani- 
cal damage,  and  moisture  concentrations. 

(4)  Investigation  of  automated  data  analysis  procedures  in  c (3 ) 
above  will  consider  both  planar  and  stereographic  data  displays." 
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2.  PROCESSING  AND  TESTING  OF  RESIN  AND  COMPOSITES 
2.1  New  Equipment 

A press  for  fabricating  flat  and  curved  laminated  plates  was  designed 
and  built  at  Texas  A&M  University.  This  press,  which  is  shown  in  Fig.  1, 
is  electrically  heated,  water  cooled,  and  designed  to  apply  temperatures 
and  pressures  in  excess  of  400°F  (204°C)  and  100  psi,  respectively.  Glass 
laminates  have  been  successfully  fabricated  and,  it  is  planned  to  initiate 
the  processing  of  graphite  composites  early  in  the  next  reporting  period. 

Nine  environmental  chambers  for  conditioning  specimens  at  various 
temperatures  and  humidities  were  also  built.  These  chambers,  two  of  which 
are  shown  in  Fig.  2,  have  an  internal  glass  cavity  6 inches  in  diameter  x 
15  inches  long.  Humidity  is  maintained  by  means  of  salt  and  acid  solutions. 
An  immersion  bath  was  constructed  for  submersing  specimens  in  water  at  con- 
trolled temperatures  up  to  the  boiling  point.  It  is  planned  to  build  five 
additional  environmental  chambers  for  use  with  existing  constant-load  test- 
ing frames;  their  design  is  described  in  Section  2.2. 

Also  acquired  during  the  first  year  are  (i)  Perkin-Elmer  thermal  test 
equipment  consisting  of  a thermo-mechanical  analyzer  (TMS-2)  and  differen- 
tial scanning  calorimeter  (PSC-2)  and  (ii)  a Micromech  microwafering  saw 
for  machining  laminates. 

2 •?  Constant  Load,  Temperature  and  Humidity  Apparatus  * 

Five  existing  lever-arm  type  constant  load  mechanical  testing  frames 
are  to  be  modified  by  the  addition  of  environmental  chambers  for  controlling 
a constant  temperature  and  humidity  during  a creep  and  recovery  test.  The 
design  of  the  chambers,  humidity  generator,  and  load  frame  modifications  have 
been  detailed.  The  current  effort  on  this  project  is  being  directed  toward 
manufacturing  and  procurement  of  the  various  components  necessary  to  make  the 
facility  operational. 

The  design  go?ls  are  to  maintain  a constant  temperature  and  humidity  en- 
vironment during  a c -eep  and  recovery  test  over  a period  as  long  as  several 
weeks.  The  temperature  ra:ige  of  interest  is  23°C  (ambient)  to  100°C  and  the 
humidity  range  is  7 - %l  ►rlaMve  humidity.  Long  term  stability  of 
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High  Temperature  Press  Figure  2.  Specimen  Conditioning  Chambers 


temperature  and  humidity  was  a major  factor  in  the  choice  of  the  system 


to  generate  the  environment.  Electric  resistance  heating  tapes  con- 


trolled by  proportional  closed  loop  temperature  controllers  will  be  used. 


Constant  humidity  will  be  maintained  by  exposing  certain  saturated  and 


non-saturated  aqueous  solutions  to  an  enclosed  air  space  at  known  temper- 


ature [1]*.  Temperature  will  be  controlled  by  placement  of  heaters  through- 


out the  system.  Humidified  air  is  to  be  circulated  by  mechanical  convection 


throughout  the  system  from  a central  humidity  generator. 


Humidity:  The  concentration  of  water  vapor  in  air  can  be  expressed 
as  parts  per  mill  ion,  percent  relative  humidity,  or  percent  of  saturation 
[2]: 


where 


C e concentration  (ppm) 
p = vapor  pressure  of  water  (mm  Ha) 


P = total  pressure  of  system  (mm  Hg) 


HR  = p 


^ x 102 


where 


HR  e relative  humidity  (%) 


Pw  = partial  pressure  of  water  at  temperature  T (mm  Hg) 


p = partial  pressure  of  water  at  saturation  at  tempera- 
ws  ture  T (rmi  Hg) 


PjP'PuJ 


u rW'  rWS' 


where 


Hp  = humidity  (K  of  saturation) 


At  room  temperature  the  partial  pressure  of  water  vapor  is  small  with  res- 
pect to  the  total  pressure.  Therefore  HR  and  Hp  are  nearly  the  same.  As 


References  are  listed  in  Section  6. 


the  temperature  increases  the  partial  pressure  becomes  more  significant 
for  a given  humidity  and  the  difference  between  percent  saturation  and 
relative  humidity  becomes  significant. 

Generating  a humidity  environment:  Production  of  standard  atmospheres 
of  known  humidity  can  be  achieved  by  several  techniques  [1,  2].  Numerous 
humidity  chambers  are  available  within  the  ranqes  of  20  - 95%  HR  and 
3 - 90°C.  However, a special  system  can  be  designed  to  operated  in  a 
2 - 98%  HR  range.  A passive  system  is  desirable  for  long  term  testing 
since  it  will  be  inherently  the  most  reliable. 

Constant  humidity  air  at  a known  temperature  can  be  maintained  in  a 
static  system  by  exposing  certain  aqueous  solutions  to  an  enclosed  air 
space  [1,  2],  Tables  of  some  of  these  solutions  are  given  in  [3,  4],  and 
selected  saturated  salt  solutions  are  listed  in  Table  1.  In  this  manner 
relative  humidity  can  be  maintained  from  7 to  98%  by  varying  the  saturated 
solution  used.  Sulfuric  acid  or  glycerine  and  water  solutions  will  also 
equil ibrate  to  a constant  humidity  in  an  enclosed  air  space  at  a known  tem- 
perature [1,  3,  4].  The  Hr  can  be  varied  by  changing  the  concentration  of 
the  solution.  However,  if  water  vapor  escapes  from  the  system  the  solution 
concentration  and  equilibrium  humidity  will  change.  With  a saturated  salt 
solution  any  loss  of  vapor  would  result  in  precipitation  of  salt  from  the 
solution  thus  maintaining  a constant  solution  concentration  and  constant 
relative  humidity.  Glycerine  solutions  cannot  be  used  above  70°C  because 
the  vapor  pressure  of  glycerine  becomes  significant  at  these  temperatures. 
For  example,  at  100°C  the  vapor  pressure  of  glycerine  is  140  mm  Ha.  Since 
the  vapor  pressure  of  sulfuric  acid  at  100°C  is  less  than  1 mm  Ha, 
sulfuric  acid  solutions  can  be  used  to  obtain  very  dry  conditions  in  the 
higher  temperature  range. 

Consequently,  in  order  to  cover  a wide  range  of  temperature  (23-100°C) 
and  relative  humidity  (7-98%),  glycerine  water,  sulfuric  acid  water,  and 
saturated  salt  solutions  will  be  used. 

Description  of  apparatus:  A system  drawing.  Fig.  3,  shows  the  compo- 
nents of  the  mechanical  testing  frame  and  humidity  generator.  The  generator 
will  be  capable  of  supplying  humidity  to  5 test  frames  although  only  one  is 
shown  for  illustration  purposes. 
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The  constant  load  test  frame  has  a single  and  compound  lever  ar- 
rangement. These  levers  provide  a mechanical  advantage  of  approximately 
9 for  the  single  lever  and  30  for  the  compound  lever.  The  maximum  load 
capacity  of  the  frame  is  6,000  pounds  on  the  specimen.  A pull  rod  extends 
through  the  base  of  the  machine  and  is  attached  to  the  specimen  within  the 
environmental  chamber.  Specimens  are  placed  in  the  chamber  by  removing 
the  bottom  cover.  Humid  air  is  circulated  through  the  chamber  by  a fan  in 
the  humidity  generator.  Insulation  surrounds  the  chamber  which  is  kept  at 
constant  temperature  by  electric  resistance  heating  elements. 

The  humidity  generator  consists  of  a shallow  container  for  the  aqueous 
solution.  This  container  is  surrounded  by  an  electrically  heated  and  in- 
sulated shell.  A small  fan  circulates  air  through  a baffle  in  the  generator 
to  promote  evaporation.  The  humid  air  is  then  circulated  to  the  test  cham- 
bers through  heated  ducts. 


2.3  Dynamic  Mechanical  Response  of  Resins  and  Composites* 
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This  portion  of  the  study  on  the  effect  of  cure  cycle  parameters  is 
concerned  witn  the  direct  monitoring  of  the  viscoelastic  complex  modulus 
of  resins  and  composites  during  and  immediately  following  the  cure  cycle. 

We  have  adapted  a piezoelectric  device  called  the  "duomorph"  for  such 
use.  Described  in  this  section  is  the  duomorph,  its  analysis,  and  the 
first  application  in  which  room-temperature  curing  epoxy  was  monitored. 

It  should  be  noted  that  a torsional  pendulum  [5,  6]  and  cone  and  plate 
devices  [7]  have  been  used  to  follow  curing  of  polymers;  but  these  tech- 
niques do  not  lend  themselves  to  in  situ  monitoring  of  laminates. 

Description  of  the  gage:  The  duomorph  consists  of  two  crystalline 
piezoelectric  plates  which  are  bonded  together;  a central  metallic  plate 
is  sometimes  used  to  alter  the  response  characteristics  of  the  assembly. 

The  duomorph  is  electrically  excited  and,  through  a study  of  its  bending 
response  when  in  contact  with  a material  of  interest,  the  viscoelastic 
properties  of  this  material  are  obtained.  The  gage  was  originally  developed 
for  the  purpose  of  measuring  in  situ  complex  modulus  when  embedded  in  the 
fuel  in  a solid  propellant  rocket  [8];  for  this  application,  which  involved 
telemetry  to  trigger  and  measure  the  gage  response,  one  piezoelectric  plate 
was  employed  as  a driver  and  the  voltage  response  of  the  other  was  used  to 
obtain  the  strain  response. 

For  precise  measurement  of  viscoelastic  properties  in  the  laboratory  it 
was  found  that  it  was  best  to  use  both  piezoelectric  plates  as  drivers  and 
employ  standard  electrical-resistance  strain  gages  to  monitor  the  bending 
of  the  assembly  [9];  see  Fig.  4.  The  two  plates  are  of  opposite  polarity 
so  that  when  they  are  excited  by  a voltage  one  plate  expands  while  the  other 
contracts.  Without  external  constraint  this  deforms  a duomorph  into  exact- 
ly the  same  shape  as  would  a uniform  bending  moment;  i.e.,  a spherical  sur- 
face. If  an  alternating  voltage  is  applied  the  curvature  oscillates  from 
one  direction  to  the  other. 

Two  pairs  of  strain  gages  are  bonded  onto  the  two  outside  surfaces 
of  the  duomorph,  and  are  used  in  a Wheatstone  bridge  to  measure  the  extent 
of  bending.  This  measurement  is  made  with  the  sensor  by  itself  and  then 

^Prepared  by  Dr.  R.  A.  Schapery  and  Mr.  B.  C.  Harbert 
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Figure  4.  Duomorph  Mounted  on  a Rubber  Block 

with  the  sensor  in  intimate  contact  with  the  material  being  tested;  the 
difference  reflects  the  stiffness  of  the  material  being  tested. 

In  most  laboratory  applications  made  to-date  circular  plates  with 
in-plane  isotropy  have  been  used,  and  the  gage  was  either  placed  in  con- 
tact with  the  surface  of  a large  mass  or  embedded  in  the  mass  of  material. 
Rubber,  solid  propellant,  soils  [8]  and  asphalt  &f0  have  been  studied  on 
other  projects. 

Work  on  the  present  AFOSR  contract  has  been  concerned  so  far  with 
developing  and  applying  a technique  to  obtain  mechanical  properties  of  a 
thin  layer  of  epoxy  resin  using  the  arrangement  shown  in  Figs.  5 and  6. 
After  the  technique  has  been  developed  and  demonstrated  for  high-temper- 
ature epoxies  employed  in  advanced  composites,  it  is  planned  to  extend 
the  concept  to  laminates,  starting  with  a unidirectional  laminate;  in 
the  latter  case,  an  orthotropic  rectangular  duomorph  which  bends  primarily 
in  the  direction  normal  to  the  fibers  will  be  used.  The  bonding  systems 
and  strain  gages  now  used  with  duomorph  are  applicable  to  temperatures 
above  400°F,  and  therefore  no  major  difficulties  are  expected  in  standard 
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Figure  5.  Cross-Section  of  Duomorph  with  Resin  Layer 


Fig.  6.  Duomorph  with  Resin  Layer  on  One  Surface 
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cure  and  post-cure  cycles. 

Analysis  of  the  Duomorph:  The  plate  configuration  in  Fig.  5 had  not 
been  analyzed  in  the  previous  studies  of  the  gage.  Therefore,  this  geo- 
metry was  analyzed  in  order  to  provide  equations  for  selecting  proper 
plate  dimensions  and  interpreting  experimental  results.  Standard  assump- 
tions of  thin-plate  bending  theory  were  made.  These  include  1)  lines  per- 
pendicular to  the  neutral  plane  remain  so  during  bending,  2)  strains  are 
small,  3)  the  plate  is  of  uniform  thickness,  and  4)  the  two  materials, 
the  sensing  disk  and  the  polymer,  are  each  isotropic  and  homogeneous  in 
the  plane  of  the  disk.  The  first  assumption  implies  that  edge  effects  are 
negligible;  it  is  assumed  that  this  will  be  true  if  the  ratio  of  the  diam- 
eter to  total  thickness  of  the  plate  is  greater  than  10.  (If  there  are  any 
designs  for  which  this  ratio  is  approached,  a finite  element  analysis  will 
be  conducted  to  check  thin  plate  theory  and  modify  the  results  if  necessary.) 

The  gage  response  quantity  of  interest  is  the  ratio  er/e0>  where  er  is 
the  difference  between  the  strains  on  the  top  and  bottom  surfaces  of  the  duo- 
morph with  the  resin,  and  eo  is  the  same  quantity  but  without  the  resin.  This 
ratio  will  be  given  first  for  elastic  materials  and  then  for  viscoelastic 
materials. 

The  normal  stress-strain  equations  for  an  element  of  elastic  resin  or 
gage  material  in  plane  stress  using  a cylindrical  coordinate  system  are 

er  = aV  + 1 (°r  'voe^  0 ) 

e0  = aV  + 1 (oe  -vor)  (2) 

where  standard  notation  for  stress  and  strain  is  employed.  The  subscripts  g 
and  r will  be  used  to  denote  mechanical  properties  of  the  duomorph  gage 
and  resin,  respectively. 

Thus 

E , v = Young's  modulus  and  Poisson's  ratio  of  the  gage  (with 
9 9 a closed  electrical  circuit) 

Er*  vr  = Young's  modulus  and  Poisson's  ratio  of  the  resin. 
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Also,  V is  the  voltage  applied  to  the  duomorph  faces  and  a is  the 
piezoelectric  expansion  coefficient;  o = o for  the  resin,  and  the 
sign  of  o for  one  plate  of  the  gage  is  opposite  that  of  the  other  one. 

The  constitutive  Eqs.  (1)  and  (2)  together  with  the  remaining 
equations  from  classical  plate  theory  for  an  unsupported  disk,[ll]  yield 


e = 3a  V 
o 


fn.  1 { 

eo  3 ( 


i }fO(h3-l)S] 


l+(h  -1 )S 
T+Th-Ds 


where 


. tr/('-vr) 


h = 1 + h. 


R ' h 


This  result  is  plotted  in  Fig.  7 for  a wide  range  of  the  two  dimension- 
less parameters,  S and  h;  the  graph  will  aid  in  selecting  appropriate 
resin  and  gage  thicknesses,  depending  on  the  resin  modulus  range  of 
interest. 

In  applications  of  the  gage,  one  is  interested  in  calculating  S, 
given  values  of  er/eQ  from  tests.  For  this  purpose  we  rearrange  Eq.  (4) 
to  find  a quadratic  equation  for  $: 

(7) 

S2|R[4(h-l)(h3-l)-3(h2-l)2]}  + S^R[4(h3+h-2)-6(h2-l )]-(h-l )}  +R-1=0 


where 


R = er/eQ 


The  largest  root  provides  the  physical  value  of  S. 
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These  results  may  be  extended  to  viscoelastic  behavior  by  means  of 
the  correspondence  principle  for  vibrating  systems  [12].  Using  complex 
notation,  the  voltage  input  is  written  as 


V = | V | e1a)t  (9) 

where  | | denotes  amplitude  and  to  is  the  frequency  (radians/time).  The 
expansion  coefficient  a is  found  to  be  a complex  function  as  a result  of 
electromechanical  interactions,  and  therefore  from  Eq.  (3), 


(10) 


!e0l  = 3 1 ot 1 1 V | and  <f>g  is  the  phase  angle  between  free  gage  expansion  and 
the  voltage.  Also,  the  bending  strain  with  resin  is  written  as 


(11) 


It  is  found  that  both  $ and  er  are  non-negative  (i.e.  the  strain  lags 
the  voltage  input  ) and  consequently  we  obtain 


R = | R | e 


-io 


R is  defined  in  Eq.  (8)  and 


(12) 


(13) 


After  | R | and  e have  been  obtained  from  duomorph  tests  at  a given 
frequency,  the  complex  version  of  S,  Eq.  (5), can  be  found  from  the 
solution  of  Eq.  (7),  in  which  R in  Eq.  (12)  is  input.  The  solution  is 
written  as 

S = | S | e '*  (14) 

Experimental  studies  indicate  that  Eg  and  are  constant  over 
the  frequency  range  employed  (0.1  - 100  HZ),  and  therefore  the  complex 
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resin  modulus  E */(l-v  *)  is  (where  an  asterisk  denotes  a complex, 
viscoelastic  property) 


= !sl 

T-v. 


ISle 


05) 


★ ★ 

The  quantities  Er  and  are  the  complex  Young's  modulus  and  Poisson  s 
ratio  of  the  resin,  respectively. 

Below  the  glass-transition  temperature  of  cured  epoxy  it  is  found 
that  v * is  essentially  a constant,  say  vr.  Therefore,  it  is  helpful  to 
rewrite  Eq.  (15)  by  introducing  the  modulus 


Thus, 


|S| 


(16) 


07) 


According  to  Eq.  (16), and  if  the  frequency  dependence  of  v * is  neglected, 
|Er*|  is  the  amplitude  of  the  complex  Young's  modulus,  Er*. 

Finally,  it  should  be  noted  that  the  maximum  value  of  S for  epoxy 
is  approximately  0.05,  as  shown  in  the  next  section.  As  a result,  we 
have  found  that  the  S term  in  Eq.  (7)  can  be  neglected;  thus  a complex 
solution  is  easily  found  analytically  in  terms  of  R = | R | e~ 1 0 by  means 
of  a small  amount  of  complex  algebra.  However,  for  completeness,  we  have 
written  a computer  program  which  solves  Eq.  (7)  exactly  in  order  to  provide 
the  complex  modulus  Er  from  gage  test  data. 

Experimental  technique:  The  duomorph  used  to  study  epoxy  consisted 
of  two  1 inch  diameter  piezoelectric  ceramic  disks  0.010  inch  thick, 
bonded  together  by  means  of  a conductive  epoxy  adhesive. 

Using  the  information  in  Fig.  7 and  Ep  a 0.5  x 10^  psi  and 

E = 11.8  x 106  psi,  it  was  estimated  that  the  thickness  ratio  h_  (resin 
9 K 

thickness  to  duomorph  thickness)  should  be  approximately  3:1.  In  preparing 
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for  casting  a layer  of  epoxy  on  one  side  of  the  duomorph  a circular  dam 
of  thin  paper  was  built  up  around  the  edge  of  the  1 inch  diameter  duo- 
morph extending  0.060  inch  above  the  surface.  The  frequency  response 
of  the  duomorph  with  the  dam  was  then  determined  with  the  assembly  resting 
on  a one  inch  thick  layer  of  polyurethane  foam.  The  position  and  support 
to  be  used  when  the  epoxy  was  to  be  placed  for  cure  . Excitation  was  from 
a sine  wave  source  at  18  volts  p-p,  and  covered  a frequency  range  of  0.1 
to  100  Hz. 

Since  the  effort  at  this  point  was  mostly  concerned  with  checking 
the  feasibility  of  following  the  cure  of  an  epoxy,  a room  temperature  cure 
system  was  selected  so  that  it  would  be  easy  to  evacuate  any  entrapped 
air  from  the  mixing  process.  A suitable  system  was  found  by  mixing  Shell 
Epon  815  with  5 parts/hundred  D.E.P.  ( Diethyl  ami noproyl amine) . A small 
amount  of  these  ingredients  was  mixed  by  hand,  evacuated,  and  then  cast 
into  the  thin  space  on  top  of  the  duomorph.  The  response  of  the  duomorph 
was  obtained  periodically  as  the  epoxy  passed  through  the  cure  process. 

Following  the  completion  of  the  cure  process  the  duomorph  with  the 
cured  epoxy  was  examined  to  check  for  thickness  variation.  It  was  found 
that  a slight  meniscus  had  formed  near  the  edge  with  the  circular  edge 
being  approximately  0.005  above  the  flat  section  part.  The  total  thick- 
ness (including  the  thickness  of  the  duomorph)  was  measured  to  be  0.087 
± .005  inch,  yielding  the  total  hR  thickness  ratio  of  3.35:1.  The  meniscus 
was  next  removed  to  the  point  where  a flat  surface  was  obtained  across  the 
epoxy  layers.  The  total  thickness  following  this  operation  was  found  to 
be  0.082  ± .002  yielding  hp  = 3.1.  Figure  6 is  a photograph  showing  the 
duomorph  with  the  epoxy  layer  after  the  surface  was  smoothed  down.  The 
"effective"  thickness  ratio  of  the  epoxy  with  meniscus  was  then  determined 
from  Eq.  (4)  to  be  3.15  by  requiring  the  resin  modulus  before  and  after 
meniscus  removal  to  be  the  same; this  value  of  h^  was  used  in  all  analyses 
of  data  obtained  during  cure. 

At  the  time  of  casting  the  layer  of  epoxy  on  the  duomorph,  similar 
shapes  were  cut  from  the  epoxy  mix  to  provide  material  with  the  same  cure 
history  for  independent  property  determination  on  the  fully  cured  material. 
The  test  apparatus  used  in  making  these  tests  consisted  of  a simple  canti- 
lever beam.  The  thin  beam  cut  from  the  epoxy  was  loaded  at  the  free  end 
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by  a weight  having  a hook  to  precisely  locate  the  point  of  application. 

Beam  deflection  at  the  free  end  was  monitored  using  a microscope  mounted 
on  a precise  traversing  stage  which  provided  linear  measurement  of  dis- 
placement to  within  ± 0.001  mm. 

Experimental  results:  The  bending  strain  amplitude  le^  and  phase  an- 
gle or  (cf.  Eq.  (11))  were  obtained  at  four  frequencies  (0.1,  1,  10,  100  Hz) 
at  ten  different  times  during  the  cure  process;  the  strain  |eQ]  and  angle 
<j>g  (cf.  Eq.  (10))  were  obtained  prior  to  the  initiation  of  curing.  Equation 
(13)  was  then  used  to  derive  R and  e at  these  times  which,  together  with 
Eqs.  (7),  (12)  and  (17)  provided  values  of  the  complex  modulus  amplitude 
|Er  | and  phase  angle  4>  . The  results  cover  a greater  than  100:1  change 
in  modulus,  and  are  shown  in  Figs.  8 and  9.  It  should  be  added  that  earlier 
direct  mechanical  calibration  of  the  duomorph  [9]  yielded 

En/(l-vJ  “ 16.8  x 106  psi  (where  vn  « 0.3). 

9 y 9 c 

All  moduli  in  Fig.  8 approach  0.45  x 10  psi  at  long  times.  This  value 

is  identical  to  that  found  independently  from  the  early-time  creep  response 
of  a beam  in  flexure  tested  at  70  hours  after  mixing  of  the  polymer  (cf.  Fig. 
8).  An  additional  beam  flexure  test  made  at  360  hours  revealed  almost  negli- 
gible creep  out  to  1000  seconds,  and  provided  again  the  same  modulus.  The 
fact  that  the  70  hour  old  resin  exhibited  significant  creep  is  very  interest- 
ing. This  behavior  implies  that  the  aging  actually  continued  well  beyond  the 
age  where  the  duomorph  data  leveled  off;  therefore,  the  longer  the  resin  was 
under  load, the  more  influence  long-term  aging  has  on  mechanical  response. 

Notice  that  the  phase  angle  data  in  Fig.  9 decreases  continuously  as 
curing  proceeds.  The  behavior  exhibited  is  very  reasonable  since  <t>  = 90° 
for  a Newtonian  fluid  and  4>  = 0°  for  an  elastic  solid. 

In  view  of  the  success  of  the  duomorph  in  following  changes  in  a room 
temperature  cure  epoxy,  it  is  planned  to  next  evaluate  it  with  the  Hercules 
3502  resin,  which  is  used  in  an  advanced  composite. 


Figure  8.  Amplitude  of  Complex  Modulus  of  Resin  vs.  Age 
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Figure  9.  Phase  Angle  of  Complex  Modulus  of  Resin  vs.  Age 
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2.4  Enthalpy  Changes  in  Epoxy  Pesin  * 

The  influence  of  cure  cycle  parameters  and  age  following  cure  on 
enthalpy  and  associated  specific  heat  changes  in  epoxy  resin  is  under 
investigation.  Emphasis  to-date  has  been  on  developing  maximum  sensi- 
tivity and  reproducibility  in  the  measurements  with  the  scanning  calori- 
meter (DSC)  since  we  are  seeking  to  measure  small  changes  in  the  thermal 
properties  produced  by  various  treatments. 

Background:  This  type  of  study  has  been  carried  out  on  other 
materials  but  the  measurements  are  difficult  in  that  the  thermal  history 
must  be  carefully  specified  [13-20].  An  additional  complication  is  that 
the  effect  of  thermal  history  is  confounded  with  the  effects  of  moisture 
content  which  is  known  to  produce  large  effects. 

The  materials  now  used  as  the  matrix  for  the  advanced  graphite  com- 
posites are  selected  to  have  an  extremely  high  value  of  the  glass  transition 
temperature,  T , to  insure  good  structural  performance  at  elevated  tempera- 
tures and  humidities;since  moisture  and  other  plasticizers  reduce  T , the 
extremely  high  value  of  T^  permits  more  reduction  before  degradation  of 
properties  becomes  a more  serious  problem.  Consequently,  the  glass  tran- 
sition in  the  dry  state  is  often  well  above  200°C  (392°F). 

The  rate  of  cure  slows  at  temperatures  below  T , and  the  glass  tran- 
sition raises  with  the  extent  of  cure;  therefore  Tg  often  rises  only  to  the 
range  of  the  cure  temperature  leaving  a number  of  unreacted  groups  in  the 
sample.  This  undercured  sample  will  then  cure  further  whenever  the  temper- 
ature and  moisture  content  become  high  enough  for  some  molecular  mobility 
to  result.  This  results  in  a glass  transition  of  the  undercured  material 
which  changes  with  time  at  any  high  temperature  which  is  needed  to  measure 
it.  One  is  limited  in  the  upper  limit  of  the  cure  temperature  by  the  tem- 
perature that  causes  chemical  degradation  of  the  sample.  It  may  happen  that 
the  glass  transition  temperature  can  be  elevated  into  this  region  of  ther- 
mal degradation  so  that,  even  in  principle,  T^  is  impossible  to  measure. 


♦Prepared  by  Dr.  J.  S.  Ham 


Hercules  Inc.  has  provided  us  with  one  gallon  each  of  3502  and 
3501-6  uncured  epoxy  resins  for  study  on  this  project.  However,  these 
resins  are  proprietary  and  no  molecular  structure  is  given.  We  will 
make  comparison  with  our  studies  of  these  resins  with  similar  studies  on 
a series  of  epoxies  whose  structure  and  composition  is  known  and  where 
reactions  can  be  completed.  These  latter  resins  are  the  series  based 
upon  DGEBA  cured  with  MDA.  DGEBA  is  the  common  commercial  epoxy  and 
comes  in  a variety  of  molecular  weights  which  can  be  used  to  control  the 
molecular  weight  between  crosslinks.  The  MDA  (methylenedianil ine)  is  a 
common  curing  agent  that  is  similar  to  the  di ami  nodi phenysulfone  probably 
used  as  a curing  agent  in  the  proprietary  resins.  These  are  high  temper- 
ature epoxy  resins  with  glass  transitions  ranging  from  95°  to  180°C 
depending  upon  the  molecular  weight  of  the  epoxy  component.  Such  materials 
have  been  studied  in  the  literature  and  have  a regular  sequence  of  proper- 
ties varying  as  the  crosslink  density  varies  over  a range  of  values  C 21 3 • 
Since  the  crosslink  density  is  one  of  the  important  parameters  used  in 
the  proprietary  aerospace  resins  to  achieve  a very  high  value  of  Tg,  we 
hope  to  be  able  to  extrapolate  from  the  standard  epoxy  materials  to  the 
properties  of  the  proprietary  resins. 

Experimental  procedures  for  casting  resins:  The  samples  are  cast  with 
some  differences  in  technique  depending  upon  the  properties  of  the  starting 
materials.  We  have  been  casting  1/8  inch  thick  slabs  using  the  standard 
techniques  of  using  two  glass  plates  held  apart  by  spacers  with  a soft  rub- 
ber tube  as  a gasket.  We  use  either  mold  release  agent  on  the  glass  or  two 
mylar  sheets  to  protect  the  glass  surface  from  the  resin.  The  mold  release 
agent  may  contaminate  the  surface  of  the  resin  but  one  obtains  smooth  flat 
surfaces.  The  mylar  sheets  sometimes  wrinkle  making  the  sample  thickness 
uneven.  For  our  thermal  measurements,  either  technique  is  possible,  but 
avoiding  the  possible  contamination  of  the  mold  release  agent  makes  the 
mylar  sheet  preferable.  We  will  need  smaller  variations  in  sample  thickness 
than  we  find  using  the  mylar  sheet  when  mechanical  measurements  are  to  be 
made  on  the  same  samples  and  so  need  to  refine  the  casting  technique. 

The  low  molecular  weight  epoxy  components  are  liquids  and  the  higher 
molecular  weight  components  are  solids.  The  MDA  is  a high  melting  solid 
but  will  dissolve  slowly  at  temperatures  below  its  true  melting  temperature. 
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The  low  molecular  weight  materials  are  mixed  with  the  MDA  to  form  a 
slurry  which  is  then  poured  into  the  molds.  The  higher  molecular 
weiqht  epoxies  are  qround  and  mixed  with  the  MDA  as  two  powders  and 
this  is  then  put  into  the  mold. 

The  mold  is  held  in  a vacuum  oven  where  we  can  watch  the  ielting 
through  a window.  The  samples  are  heated  to  between  90  and  lOC-C  vhile 
under  vacuum.  At  these  temperatures  there  is  a component  given  off 
(possibly  dissolved  water,  but  the  volume  released  seems  to  be  too  large 
for  this  to  be  the  case)  producing  a foam.  The  pressure  is  increased 
whenever  the  foam  threatens  to  overflow  the  cell  so  that,  by  varying  the 
pressure,  we  expand  and  contract  the  bubbles  and  they  gradually  work 
their  way  to  the  surface  of  the  material.  The  low  molecular  weight  resins 
and  the  proprietary  materials  become  quite  fluid  at  these  temperatures  so 
the  bubbles  can  be  worked  out  quickly.  The  higher  molecular  weight  epoxy 
components  are  viscous  and  so  considerable  patience  is  required  to  get 
the  bubbles  from  these  materials. 

While  higher  temperatures  would  make  the  material  more  fluid  it 
can  also  start  the  curing  reaction  which  is  to  be  avoided  until  the  bub- 
bles are  completely  removed.  Once  the  bubbles  are  removed,  the  samples 
are  cured  at  180°C  for  two  hours. 

The  proprietary  resins  have  a curing  reaction  triggered  at  high  tem- 
peratures by  the  decomposition  of  an  amine  complex  of  BF^.  Thus,  these 
resins  have  a broad  temperature  range  in  which  they  are  liquid  before  they 
start  curing.  We  have  cured  samples  up  to  1/8  inch  thick  and  have  had 
no  difficulty  with  thermal  runaway  which  can  happen  in  the  resin  samples 
according  to  reports.  Thermal  runaway  is  a phenomenon  where  the  reaction 
produces  heat  faster  than  it  can  be  conducted  away  and  so  the  reaction  ac- 
celerates in  the  middle  of  the  sample  to  scorch  the  material  there.  This 
can  be  avoided  when  it  is  a problem  by  initially  curing  at  a low  curing 
temperature  until  enough  of  the  reaction  is  complete  so  that  the  remainder 
will  not  get  out  of  control.  ; 

Enthalpy  measurements:  Our  first  study  is  an  attempt  to  measure  the 
relaxation  of  enthalpy  with  time  at  a temperature  below  the  glass  tran- 
sition temperature  [19,  20].  To  do  this,  the  samples  are  cured  in  the 
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scanning  calorimeter  so  that  no  further  reactions  can  be  detected  with 
repeated  heating  to  some  determined  high  temperature.  The  sample  in  its 
aluminum  holder  is  then  held  for  a determined  amount  of  time  at  some  fixed 
lower  temperature  where  relaxations  in  the  glassy  material  take  place. 

This  time  is  up  to  one  month.  The  DSC  scan  is  run  at  a high  sensitivity 
and  the  sample  is  held  above  the  glass  transition  for  5 or  10  minutes  and 
then  rerun.  The  difference  in  these  two  runs  is  the  effect  of  aging  at  the 
lower  temperature  since  the  sample  is  restored  to  its  original  condition 
by  heating  above  the  glass  transition  temperature.  The  observed  ef- 
fects are  small  and  some  typical  traces  are  shown  in  Fig. 10  but  we  have 
not  been  able  to  make  a sufficient  number  of  runs  to  get  systematic  data 
since  the  DSC  must  be  set  very  precisely  to  permit  the  slope  and  intercept 
to  be  set  appropriately  without  a trial  run. 

The  sample  in  Fig. 10  was  aged  one  month  at  room  temperature  after 
holding  above  Tq  to  erase  prior  history.  The  trace  obtained  is  the  upper- 
most curve.  After  holding  this  sample  above  T^  for  5 minutes,  the  next 
curve  was  obtained.  This  curve  has  the  normal  step  shaped  curve  typical  of 
a glass  transition.  The  difference  between  these  two  curves  is  due  to  aging 
or  enthalpy  relaxations.  (The  slope  of  the  background  has  been  selected  so 
that  the  trace  stays  on  the  scale  even  with  the  change  of  specific  heat  around 
the  glass  transition;  subtract  this  slope  from  the  data  to  obtain  the  true 
variation  of  specific  heat.  A vertical  shift  was  made  so  that  the  traces 
would  not  cross  each  other;  the  true  value  of  the  specific  heat  is  the 
same  above  the  transition  region  for  the  two  runs  shown. 

Future  work:  There  are  a number  of  correlations  to  be  made  once 

we  are  able  to  measure  the  relaxation  of  the  enthalpy  with  high  precision. 
First,  we  would  like  to  compare  the  behavior  of  the  standard  materials  with 
that  of  the  proprietary  materials  used  in  the  advanced  composites.  These 
should  have  quite  similar  behavior  except  for  T . Struik  says  that  this  en- 
thalpy relaxation  occurs  at  all  temperatures  above  that  of  the  secondary 
transition  [18].  This  should  be  similar  in  these  materials. 

Second,  we  would  like  to  look  at  the  influence  of  moisture  upon  the 
relaxation  of  the  enthalpy.  It  is  known  that  moisture  can  lower  the  glass 
transition  by  many  degrees  and  it  may  enhance  the  rates  of  enthalpy 
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Figure  10.  DSC  Scan  of  Epoxy  Resin  System 


relaxation  also.  There  is  little  information  on  this  question  in 
the  literature. 

Third,  there  is  a cluster  of  questions  concerning  the  relationship 
of  the  relaxation  in  the  enthalpy  with  the  relaxation  in  residual 
stresses  and  mechanical  measurements.  There  seems  to  be  some  correlation 
but  also  some  differences  in  the  published  data  [18,  19]. 

Finally,  after  looking  at  the  various  phenomena  related  to  the  glass 
transition,  we  plan  to  study  curing  cycle  phenomena.  This  would  involve 
the  determination  of  the  extent  of  reaction  at  various  temperatures  from 
the  heat  produced  [22].  Possibly  the  enthalpy  relaxations  will  affect 
further  curing  so  that  different  temperature  profiles  will  result  in  dif- 
ferent material  properties.  Although  a more  complete  cure  might  be  ob- 
tained by  a temperature  profile  which  increases  as  the  cure  progresses, 
this  may  not  be  desirable  in  that  the  stresses  frozen  in  during  the  cool 
down  may  produce  cracking  in  the  more  completely  cured  sample. 

These  problems  involve  the  use  of  the  DSC  as  the  major  experimental 
technique  along  with  some  mechanical  measurements.  None  are  easy  ques- 
tions to  answer  experimentally  since  they  all  require  precise  control  and 
specification  of  the  sample  history  as  well  as  high  precision  and  repro- 
ducibility of  the  measurements. 

2.5  Thermal  Expansion  Behavior  of  Epoxy  Resin* 

The  influence  of  post-cure  cooling  path  on  subsequent  mechanical 
behavior  of  epoxy  resins  is  under  investigation.  Optimizing  the  post- 
cure cooling  path  to  minimize  residual  stresses  or  other  parameters  in 
graphite/epoxy  composites  is  the  ultimate  objective  in  this  study.  The 
initial  approach  to  characterizing  physical  behavior  has  been  to  measure 
the  relative  amount  of  free  volume  retained  below  Tg  as  a function  of 
cooling  rate.  According  to  current  understandings  subsequent  aging  phe- 
nomenon would  be  associated  with  a time  dependent  reduction  of  this  free 
volume.  The  volume  measurements  are  being  made  on  a Perkin-Elmer  TMS-2 
("Thermo-Mechanical  Analyzer" ) .which  provides  the  change  in  a linear  dimen- 
sion of  a virtually  unstressed  sample  with  temperature  and  time. 
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Background:  It  is  generally  observed  that  most  epoxies  experience 

some  degree  of  change  in  mechanical  properties  with  time  (the  so-called 

"aging"  phenomenon)  subsequent  to  cooling  through  their  glass  transition 

temperature  (Tg).  As  a test  of  this  behavior  the  creep/recovery  response 

at  300°F  of  Hercules  Inc.  epoxy  resin  3502  was  first  determined  for  a 

specimen  "aged"  at  ambient  temperature  for  several  months.  The  specimen 

was  subsequently  heated  at  10°F/min.  to  425°F  (T  * 375°F),  and  then  quenched 

y 

to  300°F  where  creep/recovery  behavior  was  determined  in  two  succesive  one- 
hour  tests.  The  results  are  presented  in  Fig.  11.  A dramatic  increase  in 
compliance  immediately  after  quenching  through  T^  is  noted.  The  material 
is  still  much  softer  than  in  the  original  aged  condition  even  one  hour  af- 
ter quenching  when  the  second  creep/ recovery  test  is  initiated.  Thus,  the 
transients  associated  with  aging  after  quenching  through  Tg  appear  to  be 
sufficiently  sluggish  to  be  significant  in  their  influence  on  selecting  op- 
timum post-cure  cooling  paths. 

The  simplest  explanation  for  the  observed  aging  behavior  is  that  it  is 
associated  with  a time-dependent  densification  due  to  reduction  of  free  vol- 
ume. An  alternative  explanation  is  that  possibly  high  residual  stresses 
produced  during  quenching  are  influencing  the  mechanical  properties  through 
nonlinear  mechanisms. 

Experimental  approach:  Our  experimental  approach  has  been  to  run  a 
series  of  thermal  cycles  from  348K  to  490K  (167-422°F)  using  a constant 
heating  rate  of  lOK/min  and  cooling  rates  which  varied  from  lOK/min  to 
80K/min.  The  variation  of  specimen  thickness  with  temperature  is  moni- 
tored using  the  TMS-2.  Loads  on  the  quartz  tube  are  maintained  very  low 
to  avoid  surface  indentation  of  the  small  rectangular  specimens.  The 
linear  measurements  may  be  converted  to  volume  as  the  epoxy  is  isotropic. 
Typical  results  are  shown  in  Fig.  12. 

Discussion  of  results:  The  cure  temperature  of  the  as-received  epoxy 
was  375°F.  A well  defined  Tg  of  375°F  is  observed  on  the  first  thermal 
cycle  (not  shown  in  Fig.  12).  Subsequent  thermal  cycles  to  the  same  maxi- 
mum temperature  show  no  well-defined  T^  unless  the  maximum  temperature  used 
in  previous  cycles  is  exceeded.  This  phenomenon  of  the  effective  T equal- 
ling the  previous  maximum  thermal  cycle  temperature  (and  effectively  the 
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Figure  11.  Time  Dependent  Compliance  of  Hercules  Epoxy  3502  (using  Creep/recovery  Test) 
Measured  on  a Specimen  Previously  Aged  for  Several  Months  at  Room  Temperature 
with  Two  Subsequent  Creep/recovery  Tests  Performed  Immediately  after  Heating 
above  T and  Quenching  to  300°F  in  90s.  All  compliance  data  is  for  T=300°F. 


Thermo-mechanical  hysteresis  of  Hercules  Epoxy  3502,  as  Measured  on  the  TMS-2 
Note  the  Measured  Fractional  Change  in  Thickness  has  been  converted  to  a 
Fraction  Volume  Change. 


maximum  curing  temperature)  continues  up  to  590K  (600°F),  where  the  epoxy 
experiences  significant  degradation.  The  differential  scanning  calorimeter 
was  used  to  determine  that  each  time  the  Tg  is  exceeded  additional  curing 
occurs,  thus  raising  the  Tg  to  this  new  temperature.  The  mobility  of  the 
molecules  above  Tg  is  short-lived,  as  this  new  mobility  appears  to  allow 
further  curing  to  occur. 

While  a true  (fixed)  Tg  is  not  observed  below  the  highest  cure  temper- 
ature, the  thermal  expansion  curves  do  show  considerable  nonlinearity, 
beginning  at  about  300°F.  If  one  assumes  a two  phase  solid  begins  to 
develop  at  300°F,  say  a rubbery  phase  precipitating  in  the  glassy  phase, 
the  change  in  slope  and  a rule-of-mixtures  calculation  suggest  a 0.1  - 0.2 
volume  fraction  of  rubber.  Dr.  Kibler  at  General  Dynamics  has  recently  ob- 
served changes  in  the  viscoelastic  behavior  of[+  45]graphite/epoxy  at  about 
300°F  which  may  also  be  explained  by  such  a hypothesis. 

Because  Tg  cannot  be  exceeded  except  on  the  first  thermal  cycle  in 
epoxy  (unless  each  successive  cycle  goes  to  a higher  maximum  temperature), 
the  magnitude  of  free  volume  available  to  be  "trapped"  on  quenching  is  rather 
small.  When  the  effects  of  thermal  lag  due  to  temperature  gradients  through 
the  specimen  are  taken  into  account,  the  change  in  free  volume  fraction  af- 
ter quenching  at  50K/min.  is  estimated  to  be  no  greater  than  0.1%.  This  was 
done  by  measuring  the  width  of  the  thermomechanical  hysteresis  loop  (cf. 

Fig.  12)  for  specimens  of  thickness  0.125",  0.062",  and  0.037"  and  then 

o 

extrapolating  thermal  lag  aT  vs  (thickness)  to  zero  thickness. 

Additional  experiments  in  which  the  quenchwas  interrupted  at  300°F 
to  allow  for  aging  have  indicated  that  the  difference  in  free  volume  (be- 
tween the  lowest  and  highest  cooling  rates)  is  no  greater  than  .03%;  this 
value  is  only  1/100  of  the  free  volume  at  T^  [23].  The  small  value  for 
free  volume  change  leads  us  to  believe  that  residual  stresses,  which  also 
occur  on  quenching  but  give  no  net  free  volume,  may  be  the  more  important 
cause  of  the  observed  aging  phenomenon  for  3502  epoxy. 

Future  work:  Our  work  in  the  coming  year  will  include  additional  ex- 
periments to  more  completely  characterize  the  effect  of  various  post-cure 
cooling  paths  on  the  viscoelastic  behavior  of  the  material.  It  Is  antici- 
pated that  modelling  of  this  behavior  will  allow  us  to  better  establish 
whether  annealing  residual  stresses  or  collapse  of  free  volume  is  principal- 
ly responsible  for  the  observed  behavior. 


The  effect  of  stress  level,  temperature  and  moisture  on  the  visco- 
elastic behavior  of  graphite/ epoxy  Hercules  Inc.  AS/3502  composite 
material  [±  45  ]?s  and  [90]15  is  being  investigated.  Current  tests  are 
on  specimens  produced  by  a standard  cure  cycle.  The  results  will  pro- 
vide baseline  data  for  later  studies  on  the  effects  of  cure  cycle  para- 
meters and  post  curing. 

It  is  known  that  moisture  absorption  causes  resin  to  soften,  behaving 
as  if  it  were  effectively  closer  to  its  glass  transition  temperature  (Tg). 

We  want  to  characterize  this  behavior  to  better  understand  the  mechanistic 
basis  for  the  interaction  of  absorbed  moisture  and  mechanical  behavior; 
for  example,  the  influence  of  moisture  on  creep  and  on  residual  stress 
relaxation  is  of  particular  interest.  Our  results  to-date  include  charac- 
terization of  the  dry  composite  as  a standard  to  which  the  wet  specimen  may 
be  compared,  construction  of  environmental  chambers  to  condition  wet  speci- 
mens (cf.  Section  2.1)  and  construction  of  environmental  chambers  in  which 
the  specimens  may  be  tested  without  any  significant  change  in  moisture  con- 
tent (cf.  Section  2.2). 

Experimental  approach:  The  initial  viscoelastic  characterization  work 
has  been  done  on  dry  specimens  using  creep/ recovery  tests.  Strain  gages 
have  been  used.  More  recently,  a high  temperature,  high  sensitivity  linear 
variable  differential  transformer  (LVDT)  system  was  purchased  and  special 
grips  were  built.  This  gives  us  the  capability  of  testing  at  high  tempera- 
tures (450°F  maximum)  and  in  hot,  wet  environments  where  strain  gages  are 
not  suitable.  The  last  part  of  the  dry  specimen  characterization  was  done 
using  the  new  LVDT  system.  The  results  agreed  nicely  with  the  earlier 
strain  gage  measurements,  proving  the  reliability  of  both  techniques. 

Results:  Typical  creep  compliance  data  taken  at  75°F,  175°F,  231°F 
and  261°F  is  presented  in  Fig.  13  for  the  [t  45]2$  composite.  The  compliance 
data  were  found  to  obey  a power  law  in  time;  the  exponent  is  approximately 
independent  of  temperature  up  to  231 °F.  The  data  at  different  temperatures 
can  be  superimposed  on  a single  curve  by  choosing  appropriate  horizontal 
translations,  log  aT  values.  These  are  also  presented  in  Fig.  13.  All 
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of  the  data  on  the  continuous  curve  were  taken  at  a stress  level  of  6450 
psi;  these  results  turn  out  to  be  very  close  to  low-stress  level  data 
recently  obtained  by  General  Dynamics  on  the  T300/5208  [±  45]£S  composite 
[24]. 

Several  additional  tests  were  run  at  higher  stresses  at  175°F  and 
261°F.  The  material  is  seen  in  Fig.  13  to  be  nonlinear  at  these  higher 
Stress  levels.  The  high  temperature  data  (261°F)  seem  to  have  a slightly 
higher  slope  at  all  three  stress  levels.  A subsequent  test  of  the  specimen 
at  6450  psi  after  it  was  tested  at  10,000  psi  gave  a curve  that  extrapola- 
ted the  lower  temperature  data,  maintaining  the  same  slope  (dotted  line  in 
Fig.  13). 

These  results  may  imply  that  at  sufficiently  high  temperatures  and 
stresses,  residual  stresses  are  relaxed  so  that  the  subsequent  behavior 
at  lower  stress  levels  is  modified.  Tests  have  also  been  conducted  on 
dry  unidirectional  specimens[90]1 g.  Analysis  of  the  results  is  expected 
to  aid  in  the  understanding  of  residual  stress  effects  since  they  are  a 
minimum  in  a unidirectional  composite. 

Future  work:  The  testing  of  the  wet  specimens  is  now  beginning,  using 
specimens  soaked  to  saturation  in  75%  and  98%  relative  humidity.  Addition- 
al very  long  time  tests  on  dry  specimens  are  also  to  be  made. 
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3.  ANALYSIS  OF  COMPOSITES  AND  ADHESIVES* 


Four  papers  on  the  analysis  of  time-dependent  behavior  were  completed 
during  the  reporting  period  and  are  included  in  the  Appendix.  They  help  to 
provide  an  analytical  basis  for  predicting  effects  of  post-cure  cooling  and 
the  dependence  of  mechanical  response  and  failure  on  hot,  wet  environments. 
Work  currently  underway  is  summarized  below. 

3.1  Current  Activities 

Moisture-stress  interaction  analysis:  One  major  theoretical  study  on 
moisture  effects  in  resins  and  composites  is  now  concerned  with  theories  of 
diffusion  and  interaction  of  moisture  content  and  mechanical  response. 
Specifically,  a critical  study  of  publications  on  continuum  and  equilibrium 
and  nonequilibrium  thermodynamic  theories  is  being  conducted  for  multi -phase 
media.  The  aim  is  to  incorporate  realistic  effects  of  moisture  on  the  be- 
havior of  polymers  and  polymeric  composites.  This  investigation  should  also 
provide  guidelines  for  experimental  work  and  interpretation  of  data. 

Viscoelastic  hygrothermal  stresses  in  composite  laminates:  An  analy- 
tical investigation  has  been  initiated  to  evaluate  the  hygrothermal  stresses 
which  arise  in  a mul ti -directional  laminate  accounting  for  the  viscoelastic 
behavior  of  the  material.  The  purpose  of  the  investigaion  is  to  quantify 
the  stress  levels  which  arise  within  the  individual  plies  and  determine 
their  relation  to  the  conditioning  regimes  that  are  employed  in  various  ex- 
perimental studies  on  moisture  effects  on  composites. 

All  major  factors  which  affect  the  linear  viscoelastic  behavior  of 
composite  materials  have  been  identified  and  quantified  on  the  basis  of 
the  best  available  data.  These  include  the  swelling  strains  due  to  mois- 
ture diffusion,  and  the  temperature  and  moisture  dependence  of  the  visco- 
elastic "shift  factor"  which  governs  the  creep  compliances  of  the  composite. 

The  analytical  formulation  of  the  problem  has  been  completed  for  the 
elastic  case,  and  the  viscoelastic  formulation  for  a balanced  0°/90°  lamin- 
ate is  being  finalized.  The  moisture  diffusion,  which  is  temperature  depend 
ent, affects  the  material  compliances  differently  at  each  location  across  the 
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thickness,  since  each  such  position  is  subjected  to  a different  moisture 
history.  The  formulations  account  for  the  transient,  non-uniform  moisture 
profiles  which  prevail  within  the  laminate  for  long  times. 

The  solution  of  the  viscoelastic  problem  will  require  discretization 
of  the  analytical  expressions  in  the  time  domain  and  in  the  spatial  coordi- 
nate across  the  thickness.  The  stresses  will  then  be  evaluated  by  means 
of  an  incremental,  time-marching,  numerical  computation  (cf.  paper  by  D.  R. 
Sanders  and  Dr.  W.  E.  Haisler  in  the  Appendix). 

Optimization  of  cool-down  temperatures  in  composites:  An  optimization 
scheme  is  being  developed  to  optimize  the  cool-down  temperatures  in  order 
to  reduce  the  residual  stresses  in  multi-directional  laminates  with  specific 
attention  to  0°/90°  balanced  laminates.  The  scheme  considers  the  thermo- 
rheological  ly  complex  viscoelastic  response  of  composites  and  strives  to 
develop  an  optimal  cool-down  path  that  will  increase  the  residual  strength 
of  the  predominant  plies  at  termination  of  cooling. 

Fundamental  f racture  toughness  parameters  for  viscoelastic  materials: 

Work  is  underway  to  extend  the  J-integral  theory  for  fracture  initiation  to 
time-dependent  composite  materials.  A preliminary  version  of  the  extension 
is  given  in  the  Appendix  (cf.  paper  by  R.  A.  Schapery).  It  is  shown  that 
for  a class  of  realistic  loading  histories  and  nonlinear  constitutive  equa- 
tions one  can  use  a time-dependent  J-integral  to  characterize  fracture 
initiation.  It  is  believed  that  such  a criterion  will  prove  useful  in  charac- 
terizing fracture  behavior  of  resins  and  composites  at  elevated  temperatures 
and  humidities. 

3.2  Future  Work 

Current  analytical  studies  which  are  summarized  above  will  be  con- 
tinued. The  effort  includes  supporting  the  experimental  work  described 
in  Section  2. 

Also,  the  finite  element  modelling  (cf.  paper  by  D.  R.  Sanders  and 
Dr.  W.  E.  Haisler  in  the  Appendix)  will  be  extended.  Current  plans  are  to 
allow  the  relevant  material  parameters  to  be  functions  of  temperature, 
stress,  moisture  content  and  perhaps  other  parameters.  We  will  assume 
this  functional  dependence  is  given  in  either  analytical  form  or  tabular 
form.  The  model  will  be  utilized  to  predict  viscoelastic  behavior  of 


specimens  subjected  to  various  loading  and  environmental  conditions  in 
experiments  to  be  conducted  during  the  course  of  the  present  research. 


4.  NEUTRON  RADIOGRAPHY* 
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4.1  Summary 

Through  the  past  year  the  Nuclear  Science  Center  has  developed  an 
extremely  high  quality  neutron  radiography  facility.  The  operational 
parameters  such  as  geometric  resolution,  beam  uniformity,  film  density, 
and  neutron  flux  have  been  optimized  to  a point  where  variation  in  material 
density  of  less  than  0.5%  may  be  measured.  A variety  of  laminate  sections 
with  known  moisture  contents  has  been  radiographed  in  an  attempt  to  esta- 
blish a relationship  between  film  density  changes  for  increased  moisture 
in  the  samples.  A laminate  with  known  defects  was  also  studied. 

4.2  Beam  Characterization 

The  initial  effort  was  directed  towards  producing  neutron  radiographs 
and  inspection  techniques  adequate  for  identifying  minute  defects  or  mois- 
ture content  in  laminate  sections.  For  this  purpose  the  existing  radiographic 
facility  was  modified  with  a series  of  collimators,  scrapers  and  filters  de- 
signed to  produce  a parallel  beam  of  neutrons  of  uniform  intensity  (Fig. 14). 
The  initial  and  final  beam  parameters  are  as  follows: 


Thermal  neutron  flux 

Scattered  and  fast  neutron  content 

Gamma  intensity 

Beam  uniformity 
Beam  divergence 


Initial 

1 . 5x1 06  n/ cm2/ sec 
15% 

4.0x10®  n/cm^ 
mr 

+17% 

1.5° 


Final 

>10'  n/cnr/sec 
10% 

4.0x10®  n/cm2 
mr 

+3.5% 

1.5° 


Also,  optimum  film  developing,  handling  and  inspection  procedures  necessary 
for  reproducible  data  collections  have  been  established. 


4.3  Sensitivity  Studies 

The  application  of  neutron  radiography  to  the  study  of  composite 
materials  was  directed  in  two  areas:  (1)  the  determination  of  the  sen- 
sitivity in  measuring  absorbed  water  concentration  and  (2)  the  deter- 
mination of  the  sensitivity  in  detecting  voids,  cracks  and  delaminations. 

Neat  resin  and  graphite  and  glass  fiber-reinforced  epoxy  samples  were 
used  to  investigate  sensitivity  to  absorbed  water.  Dry  samples  were  radio- 
graphed initially  in  pairs  (2  neat,2  glass, etc),  and  the  samples  were  radio- 
graphed again  after  one  sample  of  each  material  was  saturated  with  moisture 
at  98%RH.  Evaluation  of  these  radiographs  revealed  unanticipated  geometry 
effects  and  interactions  between  the  samples  which  hampered  proper  inter- 
pretation of  the  data  collected  (see  Fig. 15).  These  effects  are  caused 
by  neutrons  scattering  from  the  sample  and  adding  to  the  background 
film  density  near  the  sample.  To  alleviate  this  problem,  samples  were  sur- 
rounded by  material  with  similar  scattering  characteristics  so  the  effect 
due  to  the  scattered  neutrons  is  uniform  in  and  around  the  sample  (see  Fig. 
16.)  Also,  due  to  the  large  lag  time  between  initial  radiography  of  the 
samples  and  radiography  of  the  samples  after  saturation  it  was  decided  to 
simulate  the  presence  of  moisture  in  dry  samples  with  a hydrocarbon  in  the 
form  of  polyethylene  sheeting.  This  procedure  has  met  with  success,  with 
the  equivalent  of  1 vi%  HpO  being  easily  detected  in  a glass  laminate  section. 
Future  work  will  pick  up  at  this  point  and  establish  a relationship  between 
the  film  density  and  percent  water  absorbed. 

In  the  area  of  void,  crack  and  delamination  detection,  studies  have 
begun  of  a test  plate*  in  which  known  defects  have  been  introduced.  Figure 
17  presents  a view  of  the  defects  which  have  been  identified  at  this  point. 

The  defects  are  referred  to  only  as  additions  or  removals  (i.e.,  an  addition 
being  a positive  change  in  the  effective  sample  thickness  and  removal  being 
a negative  change);  introduction  of  a material  which  does  not  interact  with 
the  neutron  beam  would  constitute  a removal.  A more  complete  analysis 
depends  upon  further  knowledge  of  the  relationship  between  sample  thickness 
and  film  density. 


♦Furnished  by  Dr.  G.  P.  Sendeckyj  of  the  Air  Force  Flight  Dynamics 
Laboratory. 
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FIG.  17.  DEFECT  TEST  SPECIMEN 
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Dr.  Richard  A.  Schapery  (Principal  Investigator)  - Aerospace  and 
Civil  Engineering  Departments,  Director,  Mechanics  and  Materials 
Research  Center 

Dr.  Walter  L.  Bradley  - Mechanical  Engineering  Department 

Dr.  Walter  E.  Haisler  - Aerospace  Engineering  Department 

Dr.  Joe  S.  Ham  - Physics  Department 

Dr.  Kenneth  L.  Jerina  - Civil  Engineering  Department 

Dr.  John  D.  Randall  - Nuclear  Engineering  Department,  Director, 
Nuclear  Science  Center 

Dr.  Jack  Weitsman  - Civil  Engineering  Department 

In  addition  to  these  faculty  members,  Drs.  Donald  Saylak  and 
James  S.  Noel  assisted  during  the  first  half  of  the  project  period 
in  the  development  of  the  processing  laboratory.  Their  valuable 
assistance  is  gratefully  acknowledged. 

5.2  Professional  Staff 

Mechanics  and  Materials  Research  Center: 

Mr.  Carl  Fredericksen  - Electronics  Technician 
Mr.  Bob  C.  Harbert  - Research  Engineer 
Dr.  S.  B.  Raju  - Computer  Programmer 

Nuclear  Science  Center: 

Mr.  Joseph  Taft  - Research  Engineer 

5.3  Spoken  Papers  and  Lectures  at  Conferences  Related  to  Composite 
Materials: 

R.  A.  Schapery: 

"Effect  of  Moisture  on  Flaw  Growth  in  Viscoelastic  Composites,"  Gordon 
on  Composites,  Santa  Barbara,  Ca.,  Jan.  1978. 
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"Time-Dependent  Behavior  of  Composites,"  ASTM  Symposium  on  Com- 
posite Materials:  Testing  and  Design,  New  Orleans,  La.,  March, 

1978. 

"On  the  Analysis  of  Crack  Growth  in  Non-homogeneous  Viscoelastic 
Media,"  American  Mathematical  Society  Symposium  on  Fracture 
Mechanics,  New  York,  March  1978. 

"The  Graduate  Academic  Program  in  Composites  at  Texas  ASM  Univer- 
sity," ASTM  Symposium  on  Composite  Materials,  New  Orleans,  La., 

March,  1978  and  the  International  Conference  on  Composite  Materials, 
Toronto,  Canada,  April,  1978. 

"Designing  with  Composite  Materials,"  Annual  Meeting  of  the  Society 
of  Plastics  Engineers,  Washington,  D.  C.,  April  1978. 

"Fracture  Mechanics  of  Solid  Propellants,"  ONR  International  Sym- 
posium on  Fracture  Mechanics,  Washington,  D.  C.,  September  1978. 

"Review  of  the  Composite  Materials  Research  Program  at  Texas  A&M," 
Air  Force  Mechanics  of  Composite  Review  Meeting,  Dayton,  Ohio, 
November,  1978. 

"High  Speed  Crack  Propagation  in  Solid  Propellant,"  ONR  Workshop  on 
the  Detonation-to-Deflagration  Transition  in  Solid  Propellant,  Akron, 
Ohio,  December  1978. 


J.  D.  Randall: 

"Neutron  Radiographic  Facility  at  Texas  A&M",  Air  Force  Mechanics  of 
Composite  Review  Meeting,  Dayton,  Ohio,  November,  1978. 


Y.  Weitsman: 

"A  Review  on  Moisture  Effects  in  Composites,"  Army  Materials  and 
Mechanics  Research  Center,  Boston,  Mass.,  July,  1978. 

"Mechanical  Behavior  of  Polymeric  Composite  Materials  Under  the 
Influence  of  Temperature  and  Moisture  Fluctuations,"  NSF  Workshop  on 
General  Constitutive  Relationships  for  Wood  and  Wood  Based  Materials. 
Syracuse  University,  NY,  July,  1978. 

"Stresses  in  Adhesive  Layers  due  to  Sorption  of  Moisture,"  Fourth  Con 
ference  on  Fibrous  Composites  in  Structural  Design,  San  Diego,  Cal. 
November,  1978. 
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RESIDUAL  THERMAL  STRESSES  DUE  TO  COOL-DOWN 
OF  EPOXY-RESIN  COMPOSITES 

BY 

Y.  WEITSMAN* 

Abstract 


E 


This  paper  concerns  the  residual  thermal  stresses  that  develop 
within  the  resin  of  a fiber-reinforced  composite-laminate  as  the 
material  is  cooled  from  cure  temperature  down  to  room  temperature. 

The  calculations  presented  herein  consider  the  v iscoelastic 
response  of  the  resin  and  account  for  the  temperature:  dependence 
and  the  stress  sensitivity  of  the  creep  compliance.  Comparisons  with 
linear  elasticity  indicate  that  viscoelastic  relaxation  may  reduce  the 
residual  stresses  by  about  20 %. 
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1 . Introduction 

Linear  elastic  stress  analysis  predicts  that  the  thermal  stresses, 
which  develop  within  a composite  material  due  to  typical  temperature 
excursions  experienced  by  the  structure,  may  exhaust  the  strength  of 
various  laminae  even  before  the  application  of  external  loads.  Effects 
which  were  attributed  to  high  residual  thermal  stresses  in  composite 
laminates  were  reported  by  several  investigators  [1]*  [2].  A similar 
prediction  can  be  made  for  the  micromechanical  thermal  stresses  within 
the  epoxy  resin  when  its  response  is  assumed  to  follow  linear  thermo- 
elasticity. The  assumption  of  linear-elastic  behavior  tends  to  over- 
estimate the  residual  thermal  stresses  and,  by  ignoring  the  time, 
temperature  and  stress  dependencies  of  the  resin's  response  - rules 
out  any  advantage  that  may  accrue  from  the  consideration  of  special 
temperature  histories. 

A major  portion  of  the  thermal  stress  is  introduced  into  the 
composite  laminate  during  cooling  from  cure  temperature  down  to 
room  temperature.  This  cool-down  stage  may  be  viewed  as  part  of  the 
manufacturing  process  and  therefore  lends  itself  to  careful  control. 

It  has  already  been  shown  that  the  selection  of  an  optimal  cool-down 
path  could  reduce  residual  termal  stresses  in  plexiglass  [3]  and  the 
same  idea  can  be  carried  over  to  epoxy  and  epoxy-based  composites. 

This  paper  presents  calculations  of  the  "michromechanical" 
residual  stresses  within  the  epoxy  due  to  cool-down.  Our  purpose  is 


* Numbers  in  brackets  indicate  references  listed  at  end  of  paper. 
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to  provide  a quantitative  assessment  of  the  role  of  inelastic  res- 
ponse. The  calculations  are  based  upon  a thermorheologically-complex 
viscoelastic  model,  which  includes  some  aspects  of  the  non-linearity 
that  is  attributed  to  stress  effects  on  relaxation. 

The  findings  of  this  paper  are  preliminary  since  they  are  based 
upon  incomplete  data.  Consequently  it  was  necessary  to  extrapolate 
existing  data  to  higher  temperature  ranges  and  to  infer  the  non-linear 
stress-effects  by  analogy  from  similar  materials  (glass-reinforced 
epoxy).  Perhaps  most  significantly,  due  to  the  complete  lack  of  data, 
it  was  necessary  to  omit  any  aging  effects,  which  may  be  very  prominent 
during  the  cool-down  of  a newly  created  laminate. 

Several  experimental  programs,  currently  in  progress,  should 
provide  information  which  would  enable  us  to  reassess  and  refine  the 
model  which  is  used  in  the  present  work. 

2.  Formulation 

Consider  a fiber-reinforced,  multi-directional , composite 
laminate  - h < z < h,  - « < x,y  < ~,  traction  free  at  the  boundaries 
z * +h  and  subjected  to  spatially  uniform  temperature  fluctuations. 

For  sufficiently  slow  fluctuation  of  the  ambient  temperature  and  suf- 
ficiently thin  layers  it  is  possible  to  discard  the  transient  states 
and  consider  the  quasi-static  case  of  uniform  temperature 
o(z,t)  = 0(ih,t)  = o(t).* 


♦In  typical  laminates,  the  transient  effects  may  be  Ignored  up 
to  temperature  fluctuations  of  5°K/m1n. 
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Since  the  present  analysis  aims  at  investigating  the  thermo- 
viscoelastic response  of  the  resin  let  us  assume  that  the  multi- 
directional fibers,  all  of  which  lie  in  planes z =const.,  affect 
the  behavior  of  the  resin  through  the  following  geometric  constraint 

Gx  = ey  = “f  A0  (1) 

In  (1)  df  is  the  coefficient  of  thermal  expansion  of  the  fibers, 

AO  is  the  difference  between  the  initial  and  final  temperatures, 

AO  = 0 -0,,  and  e , e are  the  in-plane  strains, 
o T ^ y 

If  we  assume  in-plane  isotropy  of  the  layer,  and  consider  linear- 
elastic  response,  we  obtain  the  following  results  [3] 


where 


a =o  = o = -cGAO 
x y 

o2  = 0 and  all  shears  vanish. 


c 


(2) 


(2*) 


In  (2)  and  (2')  G,  v are  the  shear  modulus  and  Poisson's  ratio 
of  the  resin,  and  am  is  its  coefficient  of  thermal  expansion. 

If  the  resin  responds  viscoelastically  equation  (2)  takes  the 

form  t 

o(tf)  = -c I G[s(tf)-e(t')]  6(t')dt'  (3) 

o 


where  c is  a "reduced  time" .given  by 
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u)  */  ^ 


and  G(t)  is  the  shear-relaxation  function  for  the  resin. 


(4) 
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The  purpose  of  the  subsequent  investigation  is  to  consider  and 
evaluate  the  significance  of  various  factors  which  affect  the  functlpn 
ay,  and  Investigate  the  influence  of  the  history  o(t)  on  o(ty). 

3.  A linear  Thermorheologically  Complex  Model  for  Epoxy 

On  the  basis  of  uniaxial  tension  data  obtained  for  "shell  58-68R" 
epoxy  resin  [4],  for  temperatures  ranging  between  20°F  and  160°F,  it 
is  possible  to  express  the  creep  compliance  as  follows 


D(0,t)  = Dq(g)  + Dj 


t 

<*T(0) 


q 


(5) 


where 


Do(0)  = [0.2731  +0.000574  ( O - 297)]  10"6 

Dj  = 0.01  x 10'6,  oT(0)  « exp^^  - 21 .82^  ,q  = 0.19 

In  (5),  D is  in  (kPa)"1  , the  time  t is  in  minutes,  and  the  tempera- 
ture 0 is  in  degrees  Kelvin.  Furthermore,  the  Poisson's  ratio  for  epoxy 
Is  approximately  constant  at  v = 0.35.  For  temperatures  below  glass  tran 
sltion  the  coefficients  of  thermal  expansion  are  such  that 
a = am  - cy  = 5 x 10"5m/rn/l°K. 

Consider  a composite  cooled-down  between  a cure-temperature 
©0  = 438°K  and  room-temperature  ©y  s 298°K.  For  lack  of  more  complete 
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data  assume  that  (5)  can  be  extrapolated  up  to  the  elevated  cure- 
temperature.  In  terms  of  the  imposed  thermal  strain  (am  - a^)(60  * ©( t) ) 
and  the  creep  compliance  (5)  we  obtain  the  following  integral  equation 
for  the  stress  o: 


Do(0)o(t)  + di 


/[• 


-[<1 


c(t)-c(f) 


o(t’ )dt' 


(6) 


where  a * a - a,. 

m t 

If  the  temperature  drops  linearly  between  ©Q  and  0^  over  a time  in- 
terval tf  then 

0(t)  = 0Q+Rt  with  R = (0f0Q)/t^ 
and  equation  (6)  reads 


- Rt  = Do(0q  + Rt)o(t) 

* D/[ / “t'®1577  ’/ 


sTT 


o(t’)dt' 


(7) 


where  in  (7)  aT(0(s))  = exp  - 21.82^ 


The  solution  of  (7)  was  obtained  numerically.  Following  the  scheme 
developed  in  reference  [3],  the  time  interval  t^  was  divided  into  sub- 
intervals At  = t^/n  and  each  intermediate  time  was  denoted  by  tj  = jAt. 
Denoting  ej  = 0(tj)  = ©0  + Rt^,  °j  = °(tj)  and  computing  the  n values 
of 

£j  * 5(V 

o 

and  n values  of  DQj  = DQ(oj)  = jo.2371  - 0.000574(0..  - 297)JlO"6 


the  following  expression  was  obtained  for  Oj 

°i  ■ c0j*>iUj-Ej.7),r  )•  Rtj  * 1 DlUr5M,q“M 

j-1 

'1D1H«-2>S  <VVi> 


(8) 


In  (8),  H(j-2)  is  a unit  step  function  which  vanishes  for  j < 2.  In  ad- 
dition, it  is  necessary  to  employ  the  initial  conditions  oQ  = 0,  cQ  = 0. 

Calculations  were  performed  on  an  Amdahl  470  computer.  The 
residual  stress  at  the  termination  of  cool-down  is  shown  in  Fig.  1 as 
a function  of  the  time  t^  allowed  for  the  cooling.  The  reduction  in  the 
residual  stress  with  the  increase  in  cooling-time  tf  is  due  to  visco- 
elastic relaxation.  The  value  of  o = 39430  kPa  at  tf  = 0 corresponds  to 
the  linearly-elastic  case. 


Additional  relaxation  of  stresses  occurs,  at  o = 0^.,  at  times 
t > tf.  The  calculation  of  these  stresses  requires  a slight  modifica- 
tion of  expressions  (7)  and  (8),  which  is  omitted  here.  The  additional 
relaxation,  which  is  drawn  by  the  dashed  lines  in  Fig.  1,  further  reduces 
the  residual  stresses  by  about  3 - 4%. 

4.  Optimal  Cool-Down  Path  for  the  Thermorheologically  Simple  Model 

It  has  been  shown  previously  [3]  that,  for  a prescribed  time  t^  for 
cool-down,  there  exists  an  optimal  cool-down  path  0 = o(t)  which  minimizes 
the  residual  stress  of.  A remarkable  feature  of  this  optimal  path  is  that 
it  contains  discontinuities  at  times  t = 0 and  t * t^.  It  Is  Interesting 
to  note  that  the  discontinuous  character  of  the  solution  Is  apparently 
quite  common,  since  a similar  result  was  obtained  recently  for  optimal 
strain  paths  [5].  The  sharp  discontinuity  represents  an  impractical  mathe- 
matical Idealization  and  Indicates  that  the  best  results  can  be  realized  by 
rapid  cooling  during  the  initial  and  final  stages. 


In  order  to  reduce  the  computational  effort  consider  the  approximate 
form  of  the  creep  compliance 

D(t)  = b(l+atq) 

where  b = 0.3126  x 10'6(kPa)_1  , a = 0.032  and  q = 0.19  (9) 

Equation  (9)  was  derived  from  (5),  with  0 = 75°F  (297°K)  and  when  Dq(0) 
is  approximated  by  a constant  value  DQ  = b. 

Furthermore,  it  is  more  convenient  to  employ  the  relaxation  modulus 
E(t)  instead  of  the  compliance  D(t).  This  task  is  accomplished  with 
the  aid  of  Laplace  transforms.  The  Laplace  transform  of  D(t)  yields 


p 5(p)  - b Pl+a  li^l  = ~ 
L Pq  J PE( 


PE(P) 


(10) 


where  E(p)  is  the  Laplace  transform  of  E(t)  [6].  An  approximate  in- 


version, based  upon  the  expression  [7]  E(t)  <* 


1 

pD(P) 


j yields 
p=  2t 


E(t)  = 


3.2  x 10c 


l+0.0336t 


0,19 


kPa 


(ID 


The  validity  of  (11)  was  verified  by  comparison  with  a collocation 
scheme  for  the  inversion  of  (10).  The  collocation  method  will  be  dis- 
cussed In  the  following  section  on  nonlinear  effects.  The  comparison  Is 
shown  In  Fig.  2. 

Employment  of  (11)  in  (3)  yields 


o(tf)=  -ciS(tf) 


(12) 
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To  determine  the  optimal  path  e(t)  subdivide  tf  into  n equal  portions  as 
before  and  consider  a "stair-case"  drop  in  0 

0(t)  = Ak  tk  < t < tk+1  o < k < n-1 

whereby 

0(t)  = (AR  - Ak_1 ) 6 (t-tk)  1 < k < n (14) 

and 

e(t)  = (Aq  -0Q)  6 (t)  k=0 

Substitution  of  (14)  into  (13)  gives 


Equations  (12)  and  (15)  express  the  residual  stress  o(t^)  In  terms 
of  the  unknown  intermediate  temperatures,  A^  * A(tj)  j * 0,1,...  n-1. 
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For  a stationary  value  it  is  necessary  that  3S(tf)/3Aj  = 0, 
j = 0,1... n-1.  This  results  in  n simultaneous  equations  in  the  n un- 


known A.  as  follows 

J 


V8o  szo  , ||H  ,,  V1  8Zf  , ' 1 


(16) 


k = 0,  1 ,. . . n-1 


where,  in  (16)  Z * 1,  H(i-l)  is  the  unit  step-function  which  vanishes 
for  negative  arguments,  and 


1 + a 


l/«r(Aj) 


J-i 


-6480  A 


-? 


n-1 


At 


Z>t(v 


q 


The  solution  of  (15)  and  (16)  was  obtained  by  an  iteration  technique.* 

The  resulting  optimal  path  is  shown  in  Fig.  3 for  tf  = 100  min.  It  is 
worth  noting  that  the  value  of  o(t^)  is  insensitive  to  the  details  of  the 
"tail-end"  of  the  path,  i.e.,  near  t^.  Cooling  along  the  optimal  path 
results  in  a residual  stress  o(tf)  * 32900  kPa.  This  represents  an  Im- 
provement of  about  6*  over  the  linear  cool-down  path. 


We  refrained  from  using  standard  optimization  methods,  because  those 
techniques  presuppose  a smooth  functional  which  can  be  approximated  by 
quadratic  function  (see  e.g..  Ref.  [8])  while  the  present  case  revolves 
about  a discontinuous  function. 


j 
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5.  Stress  Effects  on  Relaxation  during  Cool-Down 

We  note  that  the  residual  stresses  reach  levels  which  exceed  half 
the  ultimate  stress  of  epoxy,  au  = 52,200  kPa.  Stresses  of  such  mag- 
nitude were  found  to  influence  the  creep-response  of  epoxy-based  compos- 
ites [9],  The  stress  effect  can  be  incorporated  into  the  "shift-factor" 
ay  of  Eq.  (4)  as  follows 


aj(o.o)  = exp 


|A/gr  - (B/0r)H[o-oc(0)]  [o-oc(0)]/ou(0)  a | /17x 

i ^ "Hf 


In  (17)  0R  is  a reference  absolute  temperature,  is  a "threshold 
stress"  beyond  which  stress  starts  to  affect  aT,  c»u  is  the  ultimate  stress 
of  the  resin,  H(  ) is  a unit  step-function  which  vanishes  for  negative 
arguments  and  A,  B are  material  constants. 

From  the  data  given  in  Ref.  [9],  which  were  obtained  at  a reference 
temperature  0R  = 293°K,  and  from  recent  data  which  detailed  the  dependence 
of  oy  on  o [10],  the  following  values  were  selected 

A * 6500,  B = 2000 
ou(0)  » 55200  - 190  (0-0R) 
oc(e)  - 13800  - 47.5  (0-OR). 

where  ou  and  oc  are  in  kPa. 

A significant  reduction  in  the  computational  effort  is  achieved  if 
the  relaxation  function  E(t),  as  given  in  (11),  is  expanded  In  a sum  of 
exponentials.  To  this  purpose,  the  method  of  collocation  [7]  was  employed 
and  E(t)  was  expressed  as 


12. 


22 

E ( t ) = ^ EIe~XIt  Xj  = 10'18,  10"17,  ...103  (18) 

1=1 


The  values  of  E j , which  are  not  listed  here,  were  obtained  from 
solving  the  simultaneous  equations 

22 

^ (I+Xj/XjJ^Ej  = pE(p), 

1=1  'P=XJ 

with  pE(p)  given  in  (10).  The  form  (18)  of  E(t)  is  shown  in  Fig.  2. 

Substitution  of  (17)  and  (18)  in  (3)  and  (4)  yields 
22  t 

o(t)  l\f  exp  [-XIc(t)+XIc(f  ) J ^|-r  df  (19) 

I =1  *b 

u 

where£(u)/  ^[HsT.»(sTT 

*0  1 

with  the  form  of  aj  given  in  (17). 

Consider  a linear  cool-down  from  0Q  to  enduring  a time  interval 
tf.  Denote  the  constant  cooling  rate  by  R = (0Q-0^)/t^  and  let 
c2=aR/(l-v).  In  this  case  (19)  reduces  to 


o(t)  = c2S(t) 
where 


(20) 
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Divide  the  time  Interval  tf  into  n equal  intervals  At  and  denote 
tfc  = kAt  as  before.  The  initial  condition  is  aQ  = 0 and  the  intermediate 
temperatures  are  ©k  = 0Q  - Rt^.  If  we  employ  the  approximation  that 

oT(o,0)  = aT(ak,0k)  for  tk  < t < tk+1  (21) 

then,  as  a consequence  of  expansion  (18),  it  is  possible  to  integrate  S(t) 
analytically  over  each  sub-interval  tk  < t < tk+,.  This  leads  to  a 
substantial  reduction  in  computations. 

y 

Oenote  VR  =«T(°k.Tk)  and  Wk  = e"  k.  From  (21)  it  follows  that 
k 

ck  ■ 'i'k>  sS(ti  - w “f-i 

i=l 


Denote  further 
r-1 

t = (t  -tm  l)  W , - Cl 
rk  / j m m-1  m-1  sk 

■ m=U 


and  (r=l,..k.  k-1, 

♦rk  ■ Vk  * “r-l<‘r  ‘ W 
Straightforward  manipulations  lead  to 


.n) 


and 


- “"(‘iVk 


(22) 


°k  - C2  Sk 
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For  times  beyond  termination  of  cool-down  t^  > t^,  the  temperature 
retains  a constant  value  0^  but  the  stress  continues  to  relax.  The 
computational  scheme  (22)  remains  valid  except  that  now 

aT^°k,0k^  * (t|(  > tp). 

Results  for  a constant  cooling  rate  from  0 = 438°K  down  to 

o 

0f  = 293°K  are  shown  in  Fig.  4 for  tp  = 100  minutes,  followed  by  a relaxa- 
tion at  293°K,  and  for  tp  = 300  minutes.  Note  that  a prolongation  in  the 
cool-down  time  does  not  effectively  reduce  the  residual  stresses.  Compari- 
son with  results  based  upon  the  lineal  thermorheologically  complex  model 
shows  that  when  stress  effects  are  incorporated  into  a^.,(Eq.  (17)), the 

relaxation  is  enhanced  and  the  residual  stresses  reduce  by  an  additional 

10*. 


Finally,  consider  the  case  of  a "two-step  cooling".  In  this  case 
the  temperature  is  dropped  suddenly  from  0Q  to  0^,  held  at  0 = 0^  for 
a period  t = tp  and  dropped  sharply  again  to  0 = Op  at  t = tp.  For 
such  cooling-history  Eq.  (19)  yields 

o(t^)  = c^S(tp)  * ^ Ej  (23) 

I 

In  (23)  cm  * (0o'Gb)  • cf  = (0b  "®f) 

and 

S,  -L 

i 

where  is  defined  and  evaluated  as  before. 

Following  a numerical  scheme  similar  to  (22),  the  dependence  of 
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the  residual  stress  o(t^)  on  the  intermediate  temperature  0^  has  been 
computed  for  a "hold-time"  t^  = 100  minutes.  In  the  calculations  we 
took  0O  = 438°K  and  0^  = 293°K.  The  results  are  shown  in  Fig.  5. 

Note  that  the  residual  stress  attains  a minimum  at  the  intermediate 
temperature  ©b  = 333°K.  A similar  behavior  was  observed  experimentally 
for  epoxy  coupons  reinforced  by  unidirectional  fibers.  Those  coupons 
were  cooled-down  in  two-steps  and  a qualitative  assessment  of  the 
residual  stresses  was  obtained  by  means  of  a frozen  fringe  pattern.  It 
was  observed  that  a minimal  fringe  count  was  attained  at  an  intermediate 
value  of  temperature  ©b  (0^  > 0^  > ©Q).* 

6.  Concluding  Remarks 

It  has  been  shown  that  a linear  elastic  analysis  of  residual  stresses 
overestimates  their  magnitude  by  more  than  20%.  If  account  is  taken  of  the 
viscoelastic  response  of  the  material  it  may  be  possible  to  design  a cool- 
down path  that  will  yield  a reduced  value  of  residual  stresses.  It  seems 
that  the  details  of  the  cool-down  path  influence  the  magnitude  of  the 
residual  stresses  more  noticeably  than  the  duration  of  the  cooling  period. 

The  computations  conducted  in  the  present  paper  are  based  upon  incom- 
plete data.  It  is  intended  to  modify  the  calculations  and  analysis  in 
light  of  new  data  as  they  become  available.  It  should  be  noted  that  the 
effects  of  aging  were  not  considered  in  the  present  formulation.  These 
effects  may  be  of  prime  importance  because  it  is  conceivable  that  the 
newly  created  epoxy  may  be  undergoing  rapid  aging  while  it  is  being  cooled 
down.  An  inquiry  into  this  matter  is  now  under  way. 

♦The  author  is  indebted  to  Dr.  Frank  Crossman,  of  Lockheed  Missiles 
and  Space  Co.,  Palo  Alto,  California,  for  this  information. 
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Figure  Titles 

Fig.  1.  Dependence  of  the  Residual  Stress  o on  Cooling  Time  t^.  Linear 

Thermorheologically  Complex  Model.  Linear  (constant-Rate) 

Cool-Down  from  0 = 438°K  to  G*  = 298°K.  Dashed  lines  indicate 
o t 

Further  Relaxation  of  Residual  Stressts  at  = 298°K  After 
Termination  of  Cool-Down. 

Fig.  2.  Relaxation  Modulus  E(t)  Obtained  by  Approximate  Inversion 
of  D(t)  as  Given  in  Equation  (11)  (Dashed  line)  and  by 
Exponential -Series  Determined  Through  Collocation  (solid  line). 

Fig.  3.  Optimal  Cool-Down  Path  for  the  Thermorheologically  Compl ex  Model . 

Fig.  4.  Variation  of  Residual  Stress  o with  Time  t. 

Stress  Effected  Relaxation  (Equation  17)).  Linear  (Constant-Rate) 
Cool-Down  from  0Q  = 438°K  to  0^  = 293°K.  Solid  Line  Corresponds 
to  Cooling  Time  t^  = 100  min.  and  also  Shows  Relaxation  Beyond 
Termination  of  Cool-Down.  Dashed  Line  Corresponds  to  tf  = 

300  min. 

Fig.  5.  Two-Step  Cool-Down.  Variation  of  the  Residual  Stress  a with  the 
Intermediate-Level  Temperature  0^.  0Q  = 438°K,  0^  = 293°K, 
tf  = 100  min.  Note  minimal  residual  stress  at  0b  = 333°K. 


(in  thousands  kPa) 
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Interfacial  Stresses  in  Viscoelastic  Adhesive-Layers 
due  to  Moisture  Sorption 

by 

Y.  Weitsman* 

ABSTRACT 

This  paper  concerns  the  mechanical  behavior  of  an  epoxy  adhesive 
layer  that  is  located  between  stiff  adherends  as  the  adhesive  absorbs 
moisture  from  the  ambient  environment. 

As  the  moisture  penetrates  the  layer  by  an  extremely  slow  diffusion- 
process  the  epoxy  undergoes  a simultaneous  process  of  stress-relaxation. 
Calculations  of  viscoelastic  interfacial  stresses  were  performed  for  a 
time-dependent  response  which  is  typical  of  epoxy  and  for  a layer  geometry 
as  encountered  in  practice. 

The  results  show  that  for  exposure  to  steady  ambient  humidity  the 
viscoelastic  stresses  are  smaller  than  their  elastic  counterparts.  How- 
ever, under  fluctuating  ambient  humidity  the  viscoelastic  response  may 
cause  stress  reversals,  and  thereby  failure  modes, which  are  not  predicted 
by  elasticity  theory. 


Accepted  for  publication  in  The  International  Journal  of  Solids  and 
Structures 
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1 . Introduction 

Adhesive  bonding  forms  an  attractive  method  of  joining  structural  members. 
In  analyzing  the  stresses  which  develop  within  the  bond  it  is  necessary  to 
account  for  the  fact  that  the  adhesive  materials  respond  in  a viscoelastic 
manner  under  loads  and  their  time-dependent  behavior  is  greatly  affected  by 
moisture  and  temperature. 

When  two  adherends  are  joined  together  by  a thin  adhesive  layer  the  ad- 
hesive absorbs  moisture  from  the  ambient  environment,  at  its  exposed  edges, 
which  induces  swelling  strains  into  the  layer.  Since  the  adherends  are  much 
stiffer  than  the  adhesive  they  constrain  the  adhesive-layer  against  its  ten- 
dency to  swell,  thus  causing  the  formation  of  extremely  high  stresses  within 
the  layer.  In  this  paper  attention  is  focused  on  the  interfacial-stresses 
which  arise  at  the  interfaces  between  the  adherends  and  the  adhesive. 

Since  the  moisture  penetrates  the  layer  by  an  extremely  slow  diffusion- 
process,  the  epoxy  may  undergo  noticeable  creep  and  relaxation  while  the  • 
diffusion  process  is  in  progress.  The  main  purpose  of  this  paper  is  to 
relate  the  interaction,  which  occurs  concurrently,  between  the  two  time- 
dependent  phenomena  - moisture-diffusion  and  stress-relaxation. 

The  analysis  employes  a variational  method  and  is  inherently  approximate 
in  nature.  It  is  due  to  this  approximation  that  the  edge  singularity  in  the 
stress  field  is  replaced  here  by  finite, though  large,  values.  However,  for 
the  exceedingly  thin  layers  that  are  utilized  in  practice,  the  selected  ex- 
pressions for  the  displacement  fields  should  provide  a good  approximation. 
Furthermore,  the  approximation  should  not  detract  from  the  basic  purpose  of 
this  work  which  is  to  provide  information  about  the  relative  influence  of 
the  diffusion  and  relaxation  times. 

2.  Formulation 

Consider  a semi-infinite.  Isotropic  and  elastic  adhesive  layer,  of  thick- 
ness 2a,  between  two  semi-infinite  adherend  plates  as  shown  In  Fig.  1. 

Let  X,  Y denote  Cartesian  coordinates  and  t time.  Let  e(X,t)  represent 
the  unconstrained  swelling  Induced  by  moisture  and  y,  v be  the  elastic  shear- 

modulus  and  Poisson's  ratio,  respectively.  Assume  a state  of  plane  strain. 

«• 

Introduce  the  following  non-dimensional  quantities 
x = X/a,  y = Y/a,  e(x,t)  = e(X/a,  t) 
sx  “ °x/p’  sy  = V**.  sxy  * xxy/w 


(1) 
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Also,  denote 

* 2v/(l-2v)  , f2  = 2(l-v)/(l-2v) , f3  = 2(l+v)/(l-2v) 

Consider  now  the  following  approximate  form  for  the  displacement 
fields  within  the  adhesive- layer 

U(X,Y)  = auQ(x)  + ^ a u2(x)yZ 

13  (2) 
V(X.Y)  = a v^xjy  + 3 a v3(x)y 

In  (2)  U and  V are  components  of  displacement  in  the  X and  Y directions, 
respectively,  and  uQ,  u2,  v-j,  v3  are  dimensionless,  yet  to  be  determined 
functions  of  the  dimensionless  coordinate  x.  Note  that  for  an  induced  swell- 
ing-strain of  the  form  e(x,t)  the  displacements  U and  V possess  the  required 
symmetries  in  Y. 

Form  (2)  is  the  lowest-order  approximation  which  provides  information  on 
the  interfacial  stresses,  at  the  boundary  between  adhesive  and  adherend,  in  a 
self-contained  manner.  For  exceedingly  thin  layers  the  approximation  provides 
sufficient  accuracy. 

Expressions  (2)  are  essentially  the  same  as  used  in  reference  [1]*. 

Employing  Hooke's  law  we  obtain  the  following  non-dimensional  form  for 
stress-strain  relations: 


V 

f, 

f3 

O' 

V + } -V  1 2 

sy 

s 

f2 

f3 

0 

v!  ♦ v3 

sxy 

0 

0 

0 

1 

-e 

+ v] ' y + j v3'y3 

where  primes  designate  derivatives  with  respect  to  x. 
The  rigid  adherends  force  the  boundary  conditions: 


♦Numbers  in  brackets  indicate  references  listed  at  the  end  of  this  paper. 


u(x,  1,  t)  = 0 
v(x,  1,  t)  = 0 
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The  boundary  conditions  (4)  imply  that  the  unknown  functions 
uo,  u2,  v]t  and  v3  in  (2)  are  not  independent.  As  noted  in  [1]  the  four 
above  mentioned  functions  can  be  treated  independently  by  introducing 
suitable  Lagrange  multipliers  sn  and  st-  The  sp  and  st  are  the  dimension- 
less force-conjugates  of  v(x,  1,  t)  and  u(x,  1,  t)  respectively,  and 
represent  the  interfacial  normal  and  shear  stresses. 

Using  the  Lagrange  multipliers  defined  above  in  the  principle  of 
virtual  work  as  in  [l]5  yields  the  following  expression  for  the  variation 
of  the  internal  energy  6E 

6E  = ua 2 J |/[sxK'  ♦ J y2  6U2' ) ♦ sy(«v,  t y2  «*3)  t s^lyauj 

* y«*, ' + J y3  «»3' >1  iy  - st  («u«  ♦ \ su2)  - sst(u0  + } u2) 

- sn  (6,1  * 1 5v3>  • «„(»,  * J v3>  (5) 

Similarly,  the  variation  of  the  external  work,  6W,  is  given  by 

6W  = U&Z  ^0,y^  ^6u° + ly2  6u2^  + sXy*(°‘y^y<Svi  + \ y3  5 v3 )"j  dy  (6) 

*o  L -* 

Integration-by-parts  of  (5)  and  employment  of  (3)  yield,  upon  factoring- 
out  the  now  Independent  variations  6uq,  6u2,  6v1  and  6v3 : 

f f2  uo"  ' f3e'  + st  * 0 

~T5  f2  uo"  ‘ TC  (fl  + 1 } V1 ' + I uo  ' 1 f3  e'  + 7 st  = 0 

? , (7) 

7{fl  + V uo'  ‘ T7  vl"  “ f3e  - sn  = 0 

f*fl  + uo'  ■ 3n>  V “5  f2vl  ■ f3e  ■ sn  = 0 

§See  particularly  the  section  "Variational  Formulation  of  the  Totally  Constrained 
Case  (vQ-  0)". 
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To  these  are  adjoined  the  constraint  conditions 


uo  + 7 u2 


V1  + 3 v3 


= 0 

= 0 


(8) 


Equations  (7)  and  (8)  express  the  field  equations  of  the  problem.  The  con- 
comitant from  the  integration-by-parts  of  (5)  combines  with  (6)  to  determine 
the  boundary  conditions  at  x = 0.  This  combination,  together  with  (3),  yields 


-7  uQ  + 2 v] 1 =0 

^1  5 ^3 

2 uo'  + 7^  V1  = 2 fj  e(o’  ^ 


(9) 


Note  that  the  boundary  conditions  (9)  were  obtained  after  employment  of 
(8),  which  hold  also  at  x * 0. 

In  order  to  solve  for  u„  and  v,  we  first  eliminate  s*  between  the  first 

o I t 

pair  of  (7)  and  sR  between  the  last  pair.  Next,  we  employ  (8)  to  eliminate 
u2  and  Vg.  This  procedure  yields: 


2 f2  uQ"  + (f,  + 1)  v/  - 5 u„  -f  f,e' 


1 ' '•  ’1  ' - “o  ' 7 t3* 


(10) 


7 V " <fl  + vo'  ‘ 3f2  V1  = 0 


The  solution  of  (10)  is  expressed  with  the  aid  of  the  following  second- 
order  equation 


y f2  Z2  + + l)2  - 6f22  ' y]  Z + 15  f2  = 0 

Denote  the  roots  of  (11)  by  Zj  and  Z2  and  define 


°1 

xm 


7 » 

7 ^ *3f2 
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“1  x2  ‘ °2  X1 


(ID 
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The  solution  of  (10),  which  vanishes  as  x is 


m 00 

u°  = D L J sinh  «2(x  - s)e'(s)ds  - sinh  (x-s)e'  (s)dsl 

Y Y J 


•a-i  X -aoX 

+ A-j6  + A26 


00 

^1  = Dx1x2 J [cosh  a2^x  " s)  " cosh  a](*  " S)J  e'(s)ds 


(12) 


x 
-a,x 


A 

■ x^A^e  * x^2q 


-a„X 


In  (12)  and  A2  are  arbitrary  constants  determined  by  boundary  con- 
ditions (9). 

3.  Elastic  Solution  for  Moisture-Diffusion  Under  Constant  Ambient  Humidity 

It  has  been  observed  [2],  [3]*  that  the  saturation-moisture  level 
in  "neat"  epoxy,  as  well  as  in  epoxy-based  composites,  depends  on  the  relative 
humidity  of  the  ambient  environment.  The  dilatational  strains  that  accompany 
moisture  sorption  are  about  2/3  of  the  swelling  that  would  be  anticipated  by 
straightforward"volume  additivity". 

The  penetration  of  moisture  into  several  epoxy  resins  was  shown  to  follow 
the  classical  diffusion  equations  [2]  [3]  [4]. 

For  the  one-dimensional  diffusion  process  considered  herein,  it  follows 
that  under  constant  ambient  humidity  the  moisture- induced  swelling  is  given 
by 


e(X,t) 


* Aerfc/— — \ 
\2/kt/ 


(13) 


In  (13)  A is  a constant  which  converts  moisture  content  to  dilatational 
strain  and  k Is  the  coefficient  of  moisture  diffusion. 

2 

In  terms  of  the  nonrdimensional  time  6 = kt/a  , the  non-dimensional  form 
of  (13)  reads 


e(x,s)  = Aerfc(x/2/e) 


(14) 


♦For  most  of  this  information,  the  author  is  greatly  indebted  to  Messrs.  J.  E. 
Halklas  and  E.  L.  McKague  of  General  Dynamics  Corporation,  Fort  Worth 
Division. 
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Substitution  of  (14)  into  (12)  leads  to  closed  form  expressions  [5], [6], 
e.g. 

J sinha(x  - s)e'  (s)ds  = - ^e6a  jsinhax|l  - j ^erf + a/e) 

- erf  (— — a/e))  I - i coshax  erf(-^-  + a/e)  - erf(— — a/e)l> 

2/e  /J  L L 2/e  2/e  Jf 

00 

and  a similar  expression  for  J cosha  (x  - s)  e'(s)ds 

x 

With  the  simpler,  closed- form  version  of  uQ  and  v-|  that  replaces  all  the 
integrals  in  (12)  it  is  possible  to  evaluate  analytically  the  unknowns  A^  and 
A2  therein. 

Employing  the  boundary  conditions  (9),  with  e(o,t)  = A,  the  values  of  A 
and  Ag  are  determined  from  the  following  two-by-two  algebraic  system: 


1 


L1  L2 

A1 

°(-L2x1J2  + L1X2*V 

M1  M2 

A2 

5 f3 

■ ^ A + D(-  M-^l^  + ^2X1  ) 

(15) 


where.  In  (15) 

Lm  = 7-(2aX)m 


= 2am  + 7T  x 


m 


Nn  “-Ae 


(m  = 1,  2) 


Jm  '-“m  erf(an/*  > 


Solving  for  Aj  and  ^rom  (15)  we  obtain  the  complete  solution  for  uQ 

and  v.  In  (12),  which  enables  us  to  evaluate  the  interfacial  tractions  s 
■ n 

and  st  by  using  the  first  and  third  of  Eqs.  (7). 

After  straightforward,  though  tedious,  manipulations  we  obtain 


I 


Sp(x,fj)  = D °2X1  ^*2^2 (x » ® ) ” al ^2^ i^i(*»^)j 


2 "**1*  2 “aox 

I alPlAle  ”3  a2P2A2e 


f3  e(x,g) 


St(x,B)  = - | f2D  |a22  XjJjlx.B)  - oi|ZX2J1(x,B)J 

. 2 f (a  2a  e'“>X  + „ 2a  e"°2X\  . ^3_  ' fe 

3 M’  ’ 2 2 J 


In  (16) 


0)i2(x,s)  * - 2 e ',2(’i,2  " ’1,2) 

K1>2(x.B)  = -f  e8”1-2  (»,*.,♦  *,;2) 

, ai  ox  / \ 

where  2 = e ’ erfc  (a-|  2 ^ + X/2*^J 


V = 6 


al  ,2X 


^-a-j  2*^B  + x/2  /6  j 


Also,  pi  2 = f1  + 1 " ^ (aX^l  2 

Form  (16)  Is  the  elastic  solution  to  the  present  problem.  While  this 
form  provides  a necessary  step  toward  the  generation  of  the  viscoelastic 
solution,  because  Its  Laplace  transform  Is  readily  available.  It  cannot  be 
evaluated  numerically  for  a wide  range  of  x and  t because  J^(x,$),  K^(x,e), 
Aj  and  A2  become  numerically  unstable.*  Difficulties  of  this  sort  stem  from 
the  peculiar  geometry  of  the  adhesive  layer,  namely,  the  extremely  small 
ratio  of  thickness  to  length. 


* For  x * 3 and  t 3 1000,  with  typical  values  for  k,  a,  and  v,  J,  and  K, 

22  -22) 

Involve  products  of  numbers  of  0 (10  ) with  numbers  of  0(10  1 , thus  ex- 

ceeding the  accuracy  provided  even  by  "double  precision". 
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1 i 
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0 


fi 


The  difficulty  is  overcome  by  the  employment  of  an  approximate  form  for 
erfc  (z),  as  follows  [7]: 

.2 

erfc  z = Q(t(z))e 

(17) 


where 


Q(t)  « 


£ 


and  t(z)  = (1  + pz) 


-1 


Appropriate  substitution  of  (17)  into  (16)  yields  an  alternative, 
numerically  stable,  form  for  the  interfacial  stresses.  We  obtain 


-a2x  -o.x 

1 . /..  . 5 f3  a2P2Lle  ' a1PlL2e 

* sn(x’6)  " 3 T~2 a 


f3  erfc(x/2/e) 


+ 3 D (-02X^2 


x£ 

Je‘4s 


+ o1X2P] 


,/g  . --*-1)1-  H(x  - 2a,0)e 
2/e  ' ' 1 


0a^  -a-j  X 


+ -i -)+  h(x-Za26)Q(|a2/B  - ^d)) 
2 

Sap^-apxl  |\  * t 

H(x-2a2fj)  e & 

+ h(x  - 2a^ 6 )Q ot-j » 

+ f (»lPlM2e'“lX  -»2P2Mle’“2X)} 

|«/x,  ^7  e 48  h(x-2a26)Q^|a2/s  - 


+ J f2D 


X2 

+ H(x-2o2e)ee  2 2 j-ai2x2j^le  ^ Q(oti /B  + ~) 


- h(x-2a1e)  Q(|<*j/B  - ~ |)^  + H(x-2a1e)e  1 1 j 

. B / 2m  -°1X  2 M _a2X\  | 

+ - ^ M2e  - a2  Mie  U 


* The  numerical  value  of  a^  (j  * 1,  2,  ...5)  and  p are  given  in  section  7.1.26, 
on  page  299  of  ref.  [7].  This  approximation  is  valid  for  all  z. 


_ ■■■ 


I 


In  (18) 


««>  ■ {?  13 


»«>  ■ { '1  $ 


A = Li  M2  - > B = x2*"1  Q ( °*1  ~ X1 L2Q  (“2 


In  the  limit  t -*•  °°  the  "fully  saturated"  elastic  solution  is 


l sn(x,oo)  = "f3  []  + 7 f\  (alPlL2e  1 ' a2P2Lle  2 


1 . ,¥  x _ 5 f_3  . 2.  "V  2 . "a2Xx 

ft  st(x»“)  ” ’ 3 A ^“l  L2e  " a2  Lle  ^ 


4.  Viscoelastic  Solution  for  Moisture-Diffusion  Under  Constant  Ambient 
Humidity. 


Uniaxial-tension  experiments  [8]  on  epoxy  indicate  that  its  creep-res- 
ponse can  be  described  by  a "power  law"  form  as  follows: 


D(t;T)  - D0(T) 


+ Di  [tht] 


In  (20)  T denotes  temperature  and  a (T)  is  the  "shift-factor"  function. 
Additional  data  on  polyester  [9]  and  epoxy  [10]  indicate  that  the  creep- 
response  depends  also  on  the  moisture  content  m,  namely  D = 0 (t;  T,  m). 
However,  this  dependence  is  omitted  in  the  present  analysis.  Furthermore, 
we  shall  consider  the  case  of  uniform  temperature  for  the  entire  adhesive. 


In  most  adhesives  the  Poisson's  ratio  v remains  constant  under  a wide 
range  of  conditions.  Consequently  we  shall  consider  a viscoelastic  shear- 
compliance  of  the  form 


u(t)  - d0(l  + d^) 


where  d = 2(1  + v)  D * 
0 0 


*The  non-dimensional  viscoelastic  stresses  are  given  by  s.4  = o../d 

ij  ij  0 
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Denote  x = a2/k  and  xr  = d^1^.  Let  p be  the  variable  of  the  Laplace 
transform  in  the  time  domain  and  n = px.  Finally,  let  the  ratio  between  the 
characteristic  diffusion-time  x and  the  characteristic  creep  time  xr  be 
e = Tr/x. 

Then,  the  Laplace-transform  of  (21)  is 

-j1  pB(p)  • 1 + Lis^  (22) 

do  (6n)q 

The  reciprocity  relations  between  the  transforms  of  the  creep  and 
relaxation  functions  [11]  yield 

d pE(p)  - 1/  [ 1 + r (q  +l)/(6n)qJ  = R(en)  (23) 

4*  U 

~ Upon  insertion  of  the  explicit  expressions  for  Aj  and  A 2 into  (16), 

— the  Laplace  transforms  of  the  elastic  stresses  (denoted  by  s)  become 


jPSn(x»P)"3° 


“1X2P1  e“^n  x-e  ai  a2xlP2  e"^  x-e  2 


/n  +a. 


fn  -a. 


/n  +a. 


al P1^2e 


i4 \ aipi 

“2/ 


/iT  -a. 


a2P2M,e 
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Note  that  the  singularities  at  /n  = oj  and  /n  = <*2  are  removable. 
From  the  fundamental  equations  of  linear  viscoelasticity  [12]  it 


follows  that  the  Laplace- transforms  of  the  viscoelastic  solution  to  our 
problem( denoted  by  overbars)  are 
d. 


■j-  P sn(x»P)  = R(0n)  • i Psn(x,p) 


d 1 

p 5t(x,p)  = R(Ort)  • ^ pst(x,p) 


(25) 


In  order  to  obtain  the  viscoelastic  solution  in  the  real-time  domain  it 
is  necessary  to  invert  (25).  The  inversion  was  obtained  by  the  method  of 
collocation  [13].  Typically,  expressions  (25)  were  plotted  vs.  log  n and  a 
set  on  N points  n^(i  = 1 ... .N)  was  selected  to  cover  the  ranges  where  those 
plots  show  noticeable  curvatures.  It  was  found  that  the  shapes  of  the  above- 
mentioned  plots,  as  well  as  the  locations  of  high  curvature  regions,  depended 
on  x.  Therefore,  it  was  impossible  to  employ  any  specific  set  of  points  in 
the  inversion  and  an  appropriate  selection  had  to  be  made  for  each  value  of 
x. 

The  calculations  followed  the  scheme  of  ref.  [13].  Note,  however,  that 


for  "power-law"  creep  as  given  in  (20)  we  have  lim  D(t,T)-*»,  hence  the  equi- 

_ t-Xs 

librium  stress  vanishes.  Furthermore,  in  the  present  case  the  trans- 


form variable  is  n = pt,  whereby  the  inversion  is  obtained  in  terms  of  dimen- 
sionless time  t/t  , rather  than  real  time  t. 

Calculations  were  performed  for  an  adhesive- layer  with  the  following 
selected  properties: 


Thickness  a = 0.003",  moisture  diffusivity  k = 0.2  x 10"8  in2/sec, 


Poisson's  ratio  v = 0.45,  100%  R.  H.-  to  - swelling  conversion-factor  A = 0.03, 
creep  compliances  D * 2 x 10~8  in2/lb,  0,  * 0.2  x 10~8  1n2/lb  with  q = 0.19. 


-Q  *■  A IV  I II  I I*'  » l/| 

The  corresponding  value  of  the  shift-factor  function  was  a(T)  = 1. 


-1 


With  the  above  numerical  values  we  obtain  t = 4500  sec  and  the  value 
of  e In  (22)  Is  e * 2450. 

Numerical  evaluations  were  performed  on  an  Amdahl  470  digital  computer. 
The  results  are  shown  In  Figs.  2-9.  The  heavy  solid  lines  in  Figs.  2 - 7 
show  the  dependence  of  the  non-dimensional  viscoelastic  normal  and  tangential 


interlaminar  stresses  s„  and  s+  on  time  t at  three  fixed  stations  within  the 

n t 
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layer,  x * X/a  = 0.01,  1 and  3,  when  the  free  edge  (x  = 0)  is  exposed  to  a 
constant  R.H.  * 100*  at  t = 0.  Note  that  the  time  scale  is  logarithmic. 

For  comparison  purposes,  the  elastic  stresses  - as  given  in  (18)  - are 
shown  in  thin  solid  lines  in  Figs.  2-7.  Note  that  any  station  within  the 
layer  senses  the  build-up  and  approach  of  the  moisture  in  a gradual  manner 
and  that  a steep  increase  in  stress  occurs  after  the  passage  of  a time-span, 
which  depends  on  the  characteristic  diffusion  time  and  the  location  of  the 
station.  Beyond  that  time  the  station  notices  a "fully  saturated"  state, 
which  explains  the  leveling-off  of  the  elastic  stresses  for  longer  times. 

The  viscoelastic  values,  which  are  affected  by  relaxation,  diminish  to 
zero  with  time. 


Thus  far  all  viscoelastic  results  were  obtained  from  inversion  of  (25) 
by  means  of  collocation. 

Consider  now  an  approximate  expression  for  the  relaxation  modulus 

[14] 


E(t)  = ETt)- 


sin  irn 
irn 


(26) 


where  n = n(t)  = j:  D(t). 

Employing  (26),  the  "quasi-elastic"  (viscoelastic)  solution  is  given 
by  [13] 


sn*6" (x»t)  = F(t)sn(x,t) 
sVt-e-(x,t)  = F(t)st(x,t) 


(27) 


where,  in  (27),  sn(x,t)  and  st(x,t)  and  the  elastic  expressions  given  In 
(18)  and 

F(t)  = _ ! sin 

1 + bt^  un(t) 

n(t)  * qbtq/(l  + btq) 
with  b * 0.046  and  q = 0.19. 

It  Is  Interesting  to  note  that  for  the  particular  geometric  and  physical 
values  for  the  adhesive- layer  employed  in  the  present  problem  the  quasi- 
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elastic  approach  yielded  results  which  were  almost  indistinguishable  from 
the  collocated  inversion.  This  fortuitous  circumstance  could  not  be  anticipat- 
ed a-priori  because  the  transformed  expressions, (24)  and  (25),  do  exhibit  large 
curvatures  in  the  n domain.  In  general,  the  usefulness  of  a quasi-elastic 
method  must  be  verified  in  such  circumstances.  In  the  present  case  expres- 
sion (27)  are  certainly  easier  to  evaluate  than  (25). 


Employing  (27),  values  of  s' 


(x,t)  and  s^‘e‘(x,t)  were  calculated  for 


fixed  values  of  time,  t = 1000  secs,  and  10000  secs,  along  the  interface  be- 
tween the  adhesive  and  adherent.  These  results  are  shown  in  Figs.  8 and  9 
for  0 < x < 1.6.  (For  comparison  purposes,  the  elastic  values  are  also 
shown  in  dashed  lines). 

5.  Viscoelastic  Solution  under  Fluctuating  Ambient  Humidity. 

Thus  far  the  adhesive  layer  was  considered  to  be  subjected  to  a sudden 
exposure  to  ambient- humidity  at  time  t = 0,  with  a swelling  conversion  factor 
A. 

For  fluctuating  ambient-humidity  A = A(t)  the  interlaminar  stresses  are 
given  by  a superposition  integral  as  follows 


tF111' 


where  m = n or  t. 

In  (28)  the  kernel -functions,  under  the  integral  sign,  are  given  by  (27). 
A sample  computation  was  performed  for  a sudden  exposure  to  100%  R.H. 
followed  by  drying  at  a later  time  t = t-j.  In  this  case 

A(t)  - A[H(t)  - H(t-t])3  (29) 

whereby 

Sin’*’ (x,t)  = s^e-(x,t)  - s^e,(x,t-t1 ) (m  = n or  t)  (30) 

Results  for  t1  * 10,000  sec.  at  x * 0.01,  x * 1.0  and  x = 3.C  are  shewn  by 
the  dashed  lines  In  Figs.  2-7.  Note  that  the  superposition  of  wetting  and 
subsequent  drying  is  algebraically  additive  in  the  elastic  case.  Therefore, 


» • 
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all  the  elastic  stresses, represented  by  thin  dashed  lines  in  Fig.  2-7, 
tend  to  zero  as  t-x*  without  any  further  sign  changes.*  However,  in  the 
viscoelastic  case,  the  drying  effects  are  superimposed  on  a relaxed  wet 
state  and  this  causes  an  overcompensation  - thereby  introducing  sign-reversals 
in  the  stresses.  Eventually,  as  t-*»,  the  viscoelastic  stresses  too  tend  to 
zero. 

It  is  observed  that  viscoelasticity  predicts  detrimental  effects  that  are 
caused  by  fluctuations  in  relative  humidity  - such  as  low  cycle  fatigue  and 
interfacial  tension  - that  are  absent  in  an  elastic  analysis. 

6.  Concluding  remarks 

The  present  analysis  is  based  on  several  idealizations  which  restrict  the 
validity  of  the  results.  The  exceedingly  large  stresses  which  are  predicted  by 
linear  theories  near  the  corners  X = 0,  Y = ±a  cannot  be  borne  by  real  adhesives 
and  the  response  in  those  regions  must  be  represented  by  a suitable  non-linear 
model . 

In  addition,  the  actual  viscoelastic  response  of  resins  depends  on  the 
moisture  content,  which  introduces  yet  another  non-linearity  into  the  analysis 
of  adhesive  layers. 

The  evaluation  of  those  non-linear  effects  still  awaits  an  appropriate 
characterization  at  the  present  time.  Once  available  it  would  require  solution 
by  means  of  numerical  methods. 
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FIGURE  TITLES 


1.  The  Adhesive  Layer  Between  Rigid  Adherends  Exposed  to 
Moisture  at  X = 0. 

2.  Elastic  and  Viscoelastic  Values  of  the  Non-Dimensional 
Normal  Interlaminar  Traction  sn  at  x = X/a  = 0.01  vs.  log 
t (t  in  seconds).  Heavy  Lines  - Viscoelastic,  Thin  Lines  - 
Elastic.  Solid  Lines  - Exposure  to  a Constant  Ambient  R.H., 
Dashed  Lines  - Exposure  to  Fluctuating  Ambient  R.H.  (according 
to  Eq.  (29)). 

3.  Elastic  and  Viscoelastic  Values  of  the  Non-Dimensional  Tan- 
gential Interlaminar  Traction  s.  at  x = X/a  = 0.01  vs.  log 
t (t  in  seconds).  Heavy  Lines  - Viscoelastic,  Thin  Lines  - 
Elastic.  Solid  Lines  - Exposure  to  a Constant  Ambient  R.  H., 
Dashed  Lines  - Exposure  to  Fluctuating  Ambient  R.H.  (according 
to  Eq.  (29)). 

4.  Elastic  and  Viscoelastic  Values  of  the  Non-Dimensional  Normal 
Interlaminar  Traction  s at  x = X/a  = 1 vs. log  t (t  in  seconds). 
Heavy  Lines  - Viscoelastic,  Thin  Lines  - Elastic.  Solid  Lines 

- Exposure  to  a Constant  Ambient  R.  H.,  Dashed  Lines  - Exposure 
to  Fluctuating  Ambient  R.H.  (according  to  Eq.  (29)). 

5.  Elastic  and  Viscoelastic  Values  of  the  Non-Dimensional  Tan- 
gential Interlaminar  Traction  s.  at  x = X/a  = 1 vs.  log  t (t  in 
seconds).  Heavy  Lines  - Viscoelastic,  Thin  Lines  - Elastic. 

Solid  Lines  - Exposure  to  a Constant  Ambient  R.  H.,  Dashed  Lines 

- Exposure  to  Fluctuating  Ambient  R.  H.  (according  to  Eq.  (29)). 

6.  Elastic  and  Viscoelastic  Values  of  the  Non-Dimensional  Normal 
Interlaminar  Traction  s at  x = X/a  = 3 vs.  log  t (t  in  seconds). 
Heavy  Lines-  Viscoelastic,  Thin  Lines  - Elastic.  Solid  Lines  - 
Exposure  to  a Constant  Ambient  R.H.,  Dashed  Lines  - Exposure  to 
Fluctuating  Ambient  R.  H.  (according  to  Eq.  (29)). 

7.  Elastic  and  Viscoelastic  Values  of  the  Non-Dimensional  Tan- 
gential Interlaminar  Traction  s.  at  x = X/a  * 3 vs.  log  t (t  in 
seconds).  Heavy  Lines  - ViscoeTastic,  Thin  Lines  - Elastic. 

Solid  Lines  - Exposure  to  a Constant  Ambient  R.H.,  Dashed  Lines 

- Exposure  to  Fluctuating  Ambient  R.  H.  (according  to  Eq.  (29)). 

8.  Elastic  and  Viscoelastic  Values  of  the  Non-Dimensional  Normal 
Interlaminar  Traction  s at  Times  t s 1000  sec.  and  t = 10000 
sec.  vs.  x.  Dashed  lines  - Elastic,  Solid  Lines  - Viscoelastic. 

9.  Elastic  and  Viscoelastic  Values  of  the  Non-Dimensional  Tangential 
Interlaminar  Traction  st  at  Times  t = 1000  sec.  and  t = 10000 
sec.  vs.  x.  Dashed  lines  - Elastic,  Solid  Lines  - Viscoelastic. 
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ON  THE  ANALYSIS  OF  CRACK  INITIATION  AND  GROWTH 
IN  NONHOMOGENEOUS  VISCOELASTIC  MEDIA 
R.  A.  Schapery1 

ABSTRACT.  Prediction  of  the  rate  of  crack  growth  in  viscoelastic  media, 
as  well  as  the  time  at  which  this  growth  initiates,  is  discussed.  Emphasis  Is 
placed  on  the  particular  analytical  problems  which  are  associated  with 
composite  materials, such  as  fiber-reinforced  plastics,  and  other  nonhomogeneous 
media.  Previous  work  on  the  analysis  of  crack  speed  and  initiation  of  growth 
in  linear  viscoelastic  bodies  is  first  reviewed.  Some  new  results  are  then 
obtained  for  initiation  of  crack  growth  in  linear  and  nonlinear  viscoelastic 
media.  Certain  important  unsolved  problems  related  to  fracture  of  viscoelastic 
materials  are  indicated  in  the  concluding  remarks. 

1.  INTRODUCTION  AND  SUMMARY.  Analysis  of  the  quasi-static  fracture  of 
isotropic,  homogeneous,  linear  viscoelastic  bodies  has  received  considerable 
attention,  as  may  be  seen  in  the  reviews  by  Knauss  [1,  2].  However,  very 
little  work  has  appeared  on  the  analysis  of  more  general  classes  of  viscoelas- 
tic media. 

In  the  former  case,  there  are  many  applications  involving  traction  bound- 
ary conditions  in  which  the  stress  distribution  is  Independent  of  viscoelastic 
properties;  for  this  situation  the  stresses  are  the  same  as  in  a geometrically 
identical  elastic  body  which  is  subjected  to  the  same  boundary  tractions. 

This  feature  greatly  simplifies  the  fracture  analysis, and  leads  to  relatively 
simple  equations  for  predicting  crack  initiation  and  growth  [3-3].  On  the 
other  hand,  when  stresses  depend  on  viscoelastic  properties,  the  fracture 
analysis  is  much  more  involved.  This  dependence  normally  exists  for  non- 
homogeneous bodies;  nonhomogeneity  may  be  due  to  moisture  or  temperature 
gradients  In  polymeric  materials  for  example,  or  to  the  presence  of  more  than 
one  continuum  phase,  as  In  fiber-reinforced  plastic.  An  approximate  method 
for  calculating  crack  growth  rate  in  such  materials  in  their  linear  range  of 
behavior  Is  developed  in  [6], 
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The  analysis  in  [6)  consists  of  solving  for  stress  and  dlaplaceaent  in 
the  neighborhood  of  the  crack  tip  and  then  combining  these  results  with  an 
energy  criterion  of  failure  for  the  material  at  the  crack  tip.  This  work  is 
reviewed  in  S2  - §4.  It  will  be  seen  that  the  model  employs  a stress  distribu- 
tion which  is  bounded  at  the  crack  tip.  In  contrast  to  predictions  for  linear 
elastic  media,  the  use  of  unbounded  stresses  in  viscoelastic  media  leads  to 
physically  unacceptable  results  for  fracture  criteria  [3,  7). 

The  initiation  of  growth  is  discussed  in  J5  and  16  for  linear  and  non- 
linear viscoelastic  media,  respectively.  An  approximate  method  of  analysis  is 
described,  which  represents  an  extension  of  the  J- Integral  criterion  to  visco- 
elastic materials. 

2.  LINEAR  STRESS  AND  DEFORMATION  ANALYSIS.  In  analyzing  linear  visco- 
elastic media  with  growing  cracks,  the  problem  will  be  divided  into  two  parts. 
First,  response  to  Imposed  loads,  temperatures,  etc.  is  calculated  for  the 
case  in  which  all  cracks  are  stationary.  The  second  part  consists  of  predicting 
the  additional  response  due  to  crack  growth  alone.  Addition  of  these  two  res- 
ponses yields  the  complete  solution.  The  specified  mechanical  input  for  the 
second  part  consists  of  tractions  which  are  applied  to  the  new  surfaces 
associated  with  the  crack  growth.  If,  in  the  complete  solution,  these  surfaces 
are  to  be  traction-free,  then  the  tractions  which  are  specified  in  the  eecond 
problem  are  equal  in  magnitude,  but  opposite  in  sign,  to  tractions  in  the  flret 
problem  acting  across  suitably  defined  imaginary  surfaces;  these  surfaces 
correspond,  of  course,  to  the  actual  instantaneous  crack  surfaces.  In  some 
problems  the  actual  crack  faces  will  not  be  traction  free  due  to,  for  example. 
Interfacial  rubbing  or  injection  of  a fluid  under  pressure;  an  appropriate 
adjustment  to  the  surface  tractions  specified  in  the  second  part  must  then  be 
made. 

Inasmuch  as  this  analysis  is  concerned  with  the  prediction  of  crack 
growth,  rather  than  the  total  mechanical  state  of  the  body,  we  shall  deal  with 
only  the  second  problem.  It  should  be  added  that  the  particular  manner  in 
which  the  original  problem  is  decomposed  ale's  in  conatructlng  approximate 
local  solutions,  as  discussed  in  (6). 

General  linear,  anisotropic  viscoelastic  relationships  connecting  the 
stresses  o^  (1,  j • 1,  2,  3)  and  strains  referred  to  an  orthogonal  set 
of  Cartesian  axes  (x^)  may  be  written  in  the  form 


where  (i.J.k.*  * 1,  2,  3)  are  the  linear  viscoelastic  creep  com- 
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pliances;  these  mechanical  properties  form  a completely  symmetric  positive- 
definite  fourth-order  tensor  [ 8 J . For  notatfonal  simplicity,  explicit 
dependence  on  the  coordinates  x^  is  not  shown;  however,  all  dependent  variables 
as  well  as  the  compliances  may  vary  from  point-to-point  in  the  body.  Also, 
the  summation  convention  is  used,  in  that  a repeated  index  is  to  be  summed  out 
over  its  range. 

Crack  growth  is  assumed  to  start  at  t • o;  the  lower  limit  of  the  Integral 
in  Eq.  (1)  is  shown  as  o~,  rather  than  o,  in  order  to  allow  for  a timewise 
Jump  in  the  stresses  which  may  occur  upon  initiation  of  growth. 

It  is  not  necessary  to  Include  the  expansion  or  contraction  associated 
with  temperature  or  moisture  Inasmuch  as  these  effects  are  taken  into  account 
in  the  problem  without  crack  growth.  However,  temperature  and  moisture  can  be 
expected  to  affect  the  particular  values  of  creep  compliances  which  exist  at 
each  point  in  the  body.  In  order  to  be  able  to  draw  upon  certain  simplifying 
features  of  the  Laplace  transform, it  is  assumed  that  changes  in  the  creep  com- 
pliances due  to  timewise  variations  of  temperature,  moisture,  and  other 
physical  and  chemical  changes  are  negligible  during  short  periods  of  crack 
growth  as  determined  by  the  crack  speed  and  scale  of  the  falling  material  at 
the  crack  tip  (6] . 

Let  us  next  introduce  the  Laplace  transform  (LT), 

8tdt.  (2) 


The  quantity  f is  the  LT  of  a time-dependent  function  f;  s is  the  transform 
parameter.  In  some  cases  it  is  desirable  to  use  the  following  notation  for 
the  LT, 

L {f>  5 f (3) 

and  to  use  L *{f)  for  the  inverse  of  LT.  Also,  the  Carson  transform, 

7 5 af  (4) 
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ij 


1 (fl  + i5L) 

2 'aXl  ax^* 


(7) 


where  p is  the  density  and  u,  are  the  displacements.  The  boundary  conditions 
are  all  homogeneous  with  exception  of  those  corresponding  to  the  new  surfaces 
which  result  from  crack  growth.  (We  assume  that  there  are  no  changes  in  the 
surfaces  of  the  body  other  than  those  due  to  crack  growth;  e.g.,  interaction 
of  crack  growth  and  ablation  is  excluded.) 

Consider  next  an  elastic  boundary-value  problem.  It  is  supposed  that  the 

problem  is  in  all  respects  identical  to  that  for  the  viscoelastic  body  except 

£ 

elastic  time-independent  compliances  S^^^.say,  are  used  in  place  of  the  visco- 
elastic compliances;  the  instantaneous  crack  geometry  and  the  specified  trac- 
tions acting  on  the  new  crack  faces  arc  taken  to  be  the  same  as  for  the  visco- 
elasticity problem.  Furthermore,  let  each  S®jkl  be  equal  to  the  corresponding 
Sljkl  *or  °ne  8enerlc  real  value  of  the  transform  parameter  s;  in  general, 
these  elastic  compliances  will  be  functions  of  x^.  Assume  now  that  the  elastic 
boundary-value  problem  with  specified  growing  cracks  is  solved.  Denote  the 
time-dependent  solutions  using  a superscript  e (viz.,  o^®,  u^*,  e^*).  It 
is  important  to  recognize  that  even  though  S®^  is  a function  of  a,  through 
its  association  with  this  dependence  is  suppressed  when  solving  the 

elasticity  problem  by  using  LT  theory  or  other  means.  Indeed,  the  time- 
dependent  elastic  solutions  will  depend  on  s because  of  this  association. 

Clearly,  the  LT  of  these  elastic  solutions  satisfies  the  LT  of  the 
viscoelastic  field  equations  ( 5) — ( 7 ) . If  we  invert  these  transformed  elastic 
solutions,  and  take  Into  account  the  dependence  of  S®  on  s,  the  resulting 
time-dependent  solutions  will  obviously  satisfy  all  of  the  field  equations  of 
viscoelasticity  in  the  time  plane. 

In  some  cases  (such  as  those  discussed  in  (3)),  this  inverse  will  satisfy 
all  boundary  conditions,  including  those  along  the  growing  crack  faces,  as 
well  as  meet  the  requirement  of  continuity  of  material  ahead  of  the  crack 
tips.  When  this  happens,  the  Inverse  is  the  solution  to  the  viscoelasticity 
problem;  uniqueness  is  assumed,  not  proved.  In  general,  however,  not  all 
conditions  will  be  satisfied  by  the  inverse  of  the  elastic  transforms. 

Pursuing  this  latter  point,  we  first  note  that  the  Inverse  obviously  will 
sstlsfy  all  homogeneous  boundary  conditions  on  the  original  surfaces  slnca  the 
LT  vanishes  for  all  s.  Also,  the  requirement  of  material  continuity  ahead  of 
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the  crack  tip  la  not  violated  16],  However,  the  traction  condition  on  the 
newly  formed  crack  faces  la  not  necessarily  satisfied.  Nevertheless,  In  many 
cases  most  of  the  error  in  this  condition  can  be  eliminated  through  super- 
position of  an  elementary  solution;  the  corrected  result,  which  Is  derived  In 
[6],  then  provides  the  basis  for  the  development  of  the  equation  for  predicting 
crack  speed. 

3.  THE  CRACK  TIP  MODEL.  The  model  we  employ  is  for  the  so-called  opening 
mode  of  growth,  and  Is  shown  In  Fig.  1;  It  Is  analogous  to  that  used  by 
barenblatt  [9]  and  Dugdale  [10]  for  elastic  media.  The  (x,y,z)  coordinate 


Figure  1.  Idealised  Crack  Tip  Neighborhood* 

system  In  Fig.  1 is  stationary,  with  the  y-axls  located  at  any  convenient 
position  within  the  continuum;  y • 0 defines  the  local  crack  surface,  and  the 
z-axis  Is  perpendicular  to  the  x-y  plane.  We  assume  this  surface  Is  locally 
flat  but  do  not  restrict  its  shape  far  from  the  crack  tip  relative  to  the 
scale  of  the  failure  zone,  a.  It  is  further  supposed  that  the  radius  of  cur- 
vature of  the  crack  edge  or  tip  is  large  relative  to  this  scale.  The  shape 
of  the  crack  edge  Is  otherwise  unrestricted,  so  that  it  may  be  an  "Internal", 
"through",  or  "surface"  crack. 

The  region  between  x ■ a#  and  x • a in  Fig.  1 Is  called  the  "failure 
zone"  in  [3].  No  constitutive  assumption  Is  made  for  this  material.  Instead, 
we  simply  assume  that  the  disintegrating  material  exerts  a normal  load  per 
unit  area,  Oj,  on  the  surrounding  continuum,  where  this  continuum  is  taken  to 
be  linearly  viscoelastic;  o ^ Is  called  the  "failure"  stress  distribution.  The 
failure  process  could  therefore  be  quite  arbitrary  without  invalidating  the 
results  of  our  theory.  Indeed,  the  material  In  the  failure  zone  need  not  even 
satisfy  a continuum  hypothesis.  The  layer  of  damaged  material  on  the  crack 
faces  (x  < a^)  is  supposed  to  be  so  small  that  its  effects  on  the  motion  of 

'RfprlitM  by  pwmMoii,  from  A.  A.  Schapery.  **A  theory  of  crack  initiation  and  growth  in  viacoelaetic  media.  I.  Theoreti- 
cal development*'.  International  Journal  of  Fracture.  Vol.  II,  Feb  H7J,  pp.  I4I-IS9.  copyrighted  by  Noordhoff  International 
FuMtehlng. 
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the  continuum  la  negligible. 

It  should  be  noted  that  different  names  are  used  by  various  authors  for 
the  thin  zone  of  material  which  we  have  called  "failure  zone".  Plastic  zone, 
process  zone,  and  damage  zone  are  commonly  used  labels.  Also,  a ^ la  often 
called  the  "cohesive  stress."  We  prefer  to  use  the  word  failure  in  this 
connection  because  the  material  in  this  zone  usually  la  much  thicker  than  the 
Interatomic  spacing.it  may  contain  many  cohesive  and  adhesive  microcracks  when 
In  a composite  material,  and  of  course  it  is  ruptured  during  the  passage  of  a 
crack  tip  [3].  Also,  in  composites  and  other  materials  there  may  be  a much 
larger  zone  of  damage  (e.g.,  microcracks  or  plastic  deformation)  outside  of 
this  thin  layer;  both  of  these  zones  are  taken  into  account  in  the  nonlinear 
theory  in  $6. 

We  assume  further  that  the  neighborhood  of  the  crack  tip  P in  Fig.  1 is 
in  a state  of  plane  strain.  This  condition  is  met  provided  that  (i)  the  dis- 
tance a is  small  compared  to  the  distance  B to  the  nearest  geometric  feature; 
(11)  a Is  small  compared  to  the  radius  of  curvature  of  the  crack  edge  at  P; 

and  (ill)  the  value  of  o away  from  P can  be  neglected  in  comparison  to  o 
z y 

near  P (cf.  Fig.  2).  Just  how  small  a must  be  depends  on  the  degree  of 
nonhomogeneity  [6]. 


Figure  2.  Normal  Stress  Distribution  Along  Crack  Plane* 

The  stresses  and  displacements  due  to  the  action  of  alone  may  be  taken 
from  (3]  if  the  viscoelastic  material  in  the  crack  tip  neighborhood  is  homo- 
geneous and  isotropic;  the  superscript  f is  used  to  denofd  of-l-duced  quan- 
tities. The  normal  stress  ahead  of  the  crack  tip,  ti  > 0,  Is  the  same  for 
elastic  and  viscoelastic  media  in  this  case,  and  is  given  by 


1 fll 


‘KrnMMS  ky  pnmm*-".  <i„m  « A.  ■*»»*„»  "A  Ikaoty  of  ml  lamauo*  MS  (fowik  la  rtacoalaatk  mSh.  I ThotfS 
r*tok*Ml".  lalataau*»oal  lowaal  ..I  I fact*'.  V.S  It.  tak.  Ikl*.  rr  topynthIM  ky  NooyOkoff  latataaUuaal 


I 


ON  THE  ANALYSIS  OF  CRACK  INITIATION  AND  GROWTH 


143 


where,  with  reference  to  Fig.  2, 


5 x - a,  £ = n - x . 


(9) 


The  viscoelastic  normal  displacement  of  the  upper  crack  face  relative  to  the 
fixed  plane  y - 0 is,  for  t > t^. 


f 


CvU-T)£lfdT  . 


(10) 


where 


'■  ■ / *«“■ 


)ln 


sc  + SI 

SC-  SI 


<K'  . 


(11) 


in  which  Cy  is  the  inverse  Laplace  transform  of  the  function  Cv  found  from 


-«(1  ~ 


(12) 


E and  v are,  respectively,  the  Carson  transformed  relaxation  modulus  and 
Polaaon's  ratio  obtained  from  a uniaxial  constant  strain  test.  The  time 


t^  la  the  time  at  which  the  crack  tip  first  reaches  the  point  at  x,  and  £ 


is  as  defined  in  Eq.  (9)  with  a • a(t)  replaced  by  a(t). 

The  crack  tip  stress  and  displacement  due  to  o{,  Eqs.  (8)  and  (10),  raay 
be  considered  as  approximations  for  a material  with  temperature  and  moisture 
gradients  if  spatial  variations  of  viscoelastic  properties  are  small  in  the 
region  affected  noticeably  by  o^.  Local  Inertia  effects  due  to  Oj  have  been 
neglected;  this  simplification  does  not  necessarily  preclude  the  use  of  these 
results  In  the  presence  of  large-scale  dynamic  effects,  but  it  does  limit  our 
results  to  crack  speeds  which  are  small  compared  to  characteristic  wave  speeds 
for  the  continuum. 


By  replacing  the  creep  compliance  Cv(t)  with  another  suitably  defined 


function,  the  rasulta  in  this  section  will  apply  to  cracks  In  orthotropic 
media  and  to  cracks  between  different  isotropic  or  orthotropic  materials  (l.e. 
adhesive  cracks) [11).  However,  there  are  certain  conditions  which  must  be 
met;  the  local  crack  surfaces  must  be  parallel  to  a principal  plane  of  ort.ho- 
tropy  and,  in  the  case  of  adhesive  cracks,  the  two  different  media  must  be 
incompressible. 

The  stress  and  displacement  distributions  due  to  o,,  Fqs.  (8)  and  (10), 
are  to  be  added  to  those  due  to  the  unloading  of  the  crack  faces  during 
propagation;  these  latter  solutions  are  derived  by  using  the  procedure 
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described  in  S2  for  the  case  o{  » 0.  It  is  not  necessary  to  superpose  the 
solution  for  the  body  without  cracks;  the  associated  crack  opening  displace- 
ment is  zero  and  therefore  this  solution  does  not  contribute  to  the  energy 
criterion  of  failure  used  in  $4. 

Now,  assume  o^  « 0 and  that  the  normal  stress  o^°,  say,  along  the  crack 
plane  is  independent  of  viscoelastic  properties;  in  the  neighborhood  of  the 
crack  tip  0y°  becomes  (3), 


<1  ” ° 


(13) 


where  Kj  is  the  stress  intensity  factor;  this  result  is  identical  to  that  for 
an  elastic  material.  The  corresponding  displacement  is 


t 

-i-  f 

m J 


3(K./e) 

cv(t  -t)  — at — dT>  1 > £i 


(U) 


The  stress  obtained  by  adding  Eqs.  (8)  and  (13)  is  bounded  at  the  crack 
tip.  - 0.  if 


h -Vi  / Lt  « 

•o  Jl 


(15) 


Define 


f = oJa 
o f o 


n = €/a 


(16) 


Then,  Eq.  (15)  may  be  rewritten  as 


T Kr 

. JL  — L. 

2 o I 
o o_ 


(17) 


where  a is  the  failure  stress  at  (.  = o and 
° 


•n 


dn 


(18) 


If  f is  a function  of  onlv  n (i.e.,  if  the  shape  of  the  stress  distribution 
o 

o„  is  constant)  then  I is  a constant.  If  in  addition  a is  Independent  of 
f O o 

a,  then  Eq.  (17)  provides  an  explicit  result  for  the  length  of  the  failure 
zone.  Use  of  Eq.  (16)  in  the  prediction  of  a for  more  general  cases  is  dis- 
cussed in  (5). 

Let  us  now  add  Eqs. (8)  and  (13)  for  stress,  Eqs.  (10)  and  (14)  for  dis- 
placement, and  eliminate  Kj  by  means  of  Eq.  (15).  There  results. 
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where 
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f °f 

/ — ^2  0 

•'o  ^ «+V 
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^ /Cv(t-T>s t2>  • c>  C1 


£(x,t)  s - 53 f + 2^  d*' 


(19) 


(20) 


(21) 


Let  us  now  consider  the  case  in  which  is  a function  of  viscoelastic 
properties.  The  procedure  outlined  in  52  leads  to  an  incorrect  traction  dis- 
tribution on  the  crack  faces  if  the  given  of  Is  used  in  the  bounded  elastic 
solution.  However,  we  may  instead  use  another  failure  stress  distribution 
o^e,  say,  in  the  elastic  analysis.  This  function  is  selected  so  as  to 
produce  the  proper  traction  condition  on  the  viscoelastic  crack  faces,  which 
leads  to  an  Integral  equation  for  . The  solution  of  this  problem  is 
lengthy,  and  since  it  is  detailed  in  [6]  it  will  not  oe  given  here.  This 
analysis  results  in  a rather  simple  approximate  equation  for  predicting  crack 
speed,  as  discussed  in  the  next  section. 


4.  PREDICTION  OF  CRACK  CROWTH  RATE.  Define  T such  that  2TdAf  is  the 
mechanical  energy  required  to  fall  an  element  of  material  of  Infinitesimal  area 
dA^;  T,  which  is  called  the  fracture  enetav.  does  not  include  the  energy  of 
the  element  at  the  starting  point  for  failure,  t » 0.  Thus,  at  some  fixed 
material  point,  x, 


r - 


(22) 


where  v^  is  one-half  the  crack  opening  displacement  at  the  point  of  rupture, 
C “ a (cf.  Fig.  1).  The  governing  aquation  for  crack  growth  is  derived  by 
substituting  into  Eq.  (22)  the  displacement  v,  such  as  that  in  Eq.  (20)  or 
that  in  [6, tq,  (39)]  for  a more  general  class  of  problems.  In  the  former 
situation  we  assume  that  crack  speed, a,  is  essentially  constant  during  the 
time  a/a,  which  is  the  time  Interval  in  which  Che  crack  tip  propagates  an 
amount  a;  it  is  further  assumed  that  o^  and  a are  essentially  timewise 
constant  during  the  same  period.  For  the  latter  case,  which  is  studied  in 
[6],  this  time  Interval  may  have  to  be  somewhat  larger,  possibly  as  large  as 
10  a/a.  Additional  simplification  of  the  analysis  is  achieved  by  accounting 
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for  the  fact  that  the  second  derivative  of  the  viscoelastic  creep  compliances 
and  their  logarithms  with  respect  to  the  logarithm  of  time  is  typically  quite 
SMll.  There  results,  finally  (6), 


r 


I Cv(ta> 


(23) 


where 


Also, 


t = a/ 3a  . 


(24) 


* I “<*«>  Ki'<3ta> 


1 

2 


(25) 


is  the  geometric  mean  of  so-called  quasi-elastic  stress  intensity  factors; 
the  quantity  Kj*(.)  is  the  stress  intensity  factor  due  to  crack  growth  (as 
defined  in  12)  in  a material  which  is  elastic, except  in  place  of  elastic  com- 
pliances there  are  viscoelastic  creep  compliances  that  are  evaluated  at  the 

"local”  times  t and  3t  . The  length  of  the  failure  zone  is 
a u 


a 


(26) 


When  stresses  are  independent  of  the  viscoelastic  creep  properties,  Eqs. 
(23)  and  (26)  reduce  to  those  derived  in  (4),  Thus, 


r*i  W KI2  <27> 

whers  a is  that  in  Eq.  (17). 

The  quantities  r and  o 1 are  defined  by  the  behavior  of  the  material  in 
o o 

the  failure  zone,  and  as  such  are  fracture  properties.  In  general  they  are 
functions  of  crack  speed,  although  they  appear  to  be  constant  for  the  poly- 
meric material  studied  in  |5).  Given  these  properties  and  the  viscoelastic 
creep  compliances,  Eqs.  (23)  - (26)  may  be  solved  co  predict  the  instantaneous 
crack  speed.  Experimental  determination  of  the  fracture  properties  is  dis- 
cussed in  (5]. 

5.  FRACTURE  INITIATION  IN  LINEAR  MEDIA.  The  approach  used  in  [4]  to 
predicting  initiation  of  crack  growth  consists  of  applying  Eq.  (22)  to  the 
initial  crack  tip,  which  is  located  at  x - 0;  the  initial  length  of  the 
failure  zone  is  assumed  to  be  zero.  The  Instantaneous  displacement  v at 
x • 0 and  the  failure  zone  length  are  calculated  from  Eqs.  (20)  and  (17), 
respectively.  For  the  case  in  which  is  constant  in  space  and  time  and 
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(28) 


r “ £ Ki2(ti)  cv2)(ti)*  (29) 

(2) 

unere  C^is  a so-called  secant  compliance, 

Cv2)(t)  5 ~2  /cv(t  • 

*■1  o” 

The  quantity  t^  is  the  time  at  which  the  mechanical  work  done  on  the  failure 
zone  at  x ■ 0 la  equal  to  the  fracture  energy  T.  Inasmuch  as  the  deformation 
history  of  the  failure  zone  In  the  Initiation  problem  Is  different  from  that 
discussed  previously  for  crack  growth  rate,  the  fracture  energies  In  the  two 
cases  are  not  necessarily  the  same. 

These  results  were  derived  in  (4)  for  the  case  in  which  stresses  are  In- 
dependent of  viscoelastic  properties.  However,  they  also  apply  to  the  more 
general  class  of  problems  covered  in  [6]  as  long  as  the  conditions  stated  In 
S3  are  met  and  the  stress  Intensity  factor  Is  Interpreted  as  that  due  only  to 
growth  of  a crack  by  the  amount  a.  For,  if  the  material  In  the  neighborhood 
of  the  tip  is  homogeneous,  the  stresses  due  to  tractions  on  the  crack  faces 
In  the  Interval  o<x<a  are  Independent  of  viscoelastic  properties;  of  course, 
these  surface  tractions  may  themselves  depend  on  viscoelastic  properties  since 
they  are  for  a globally  nonhomogeneous  body  without  cracks. 

The  secant  compliance,  Eq.  (30),  reduces  to  the  creep  compliance,  Cy(t), 
when  Kj  Is  constant  during  the  Interval  o < t < t^.  In  this  case 

r ■ i *i2<Vcv(ti>  (31) 

There  are  many  other  cases  In  which  Eq.  (31)  is  a good  approximation;  viz., 
c<2>  ■ C If  C is  a sufficiently  weak  function  of  time  or  If  the  second 

V V V , 

derivative  of  Cy  and  Kj  with  respect  to  the  logarithm  of  time  Is  small  (8]. 

Recall  that  o{  was  assumed  to  be  spscewise  and  timewise  constant  In 
obtaining  Eqs.  (28)  and  (30).  Whether  or  not  this  Is  a reasonably  good 
assumption  depends  on  the  particular  material  of  interest.  If,  Instead  of 
assuming  o^  is  timewise  constant,  we  assume  a is  constant,  then  It  Is  found 
Eq.  (29)  must  be  replaced  by 


(30) 


the  stress- Intensity  factor  does  not  decrease  in  time,  we  find 


and  then  obtain  an  implicit  equation  for  "initiation  time",  t^: 


A 


I 


l 
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r * i W f\(t  -*>  ar dT  <32) 

o “ 

Equation  (28)  still  applies,  but  oq  will  vary  with  time  in  order  to  satisfy 
this  equation.  The  assumption  of  a constant  value  for  a may  be  a valid 
approximation  in  some  cases,  such  as  for  composite  materials  in  which  the  size 
of  the  damage  zone  is  closely  related  to  the  scale  of  the  phases.  However, 
very  little  is  known  at  this  time  about  the  time-dependence  of  the  failure 
zone  size  in  real  materials. 

6.  FRACTURE  INITIATION  IN  NONLINEAR  MEDIA.  In  this  section  we  shall 
describe  a tentative  method  for  characterizing  and  predicting  the  initiation 
of  crack  growth  in  nonlinear  viscoelastic  materials.  The  zone  of  intense 
damage  and  rupture  is  still  assumed  to  be  the  thin  strip  of  length  a in  Fig. 

1.  However,  outside  of  this  zone  the  material  may  be  nonlinear  and  suffer 
damage.  This  approach  is  motivated,  in  part,  by  the  almllarity  of  the  linear 
viscoelastic  result  for  fracture  initiation  (cf.  Eq.  29)  and  that  for  an 
elastic  material  (C  - constant),  as  well  as  the  success  of  the  J-lntegral 
criterion  for  elastic-plastic  and  nonlinear  elastic  media  [e.g.  12-14]. 

In  order  to  illustrate  the  method,  consider  the  viscoelastic  test  speci- 
men In  Fig.  3 which  contains  an  initial  crack  of  length  a.  The  applied  force 
Is  given  as  F - Rt",  where  R and  n are  positive  constants;  we  suppose  that 
several  different  specimens  are  tested  using  different  values  of  R and  a, 
but  with  the  same  value  of  n. Figure  4 shows  displacement-time  curves  for  e 
viscoelastic  material  prior  to  crack  propagation.  From  these  results  the 


Figure  3.  Plate  Specimen  for 
Fracture  Test 


Figure  4.  Corner  Displacement 
Response  for  One  Initial  Crack 
Length 
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Isochronal  force-displacement  curves  In  Fig.  5 can  be  constructed. 


0 


Figure  5.  Force  Displacement  Curves 
Constructed  from  Figure  4. 


If  the  material  were  elastic  or  elastic-plastic  there  would  be  no  dependence 
on  time,  and  it  would  be  possible  to  use  such  data  from  specimens  with  dif- 
ferent Initial  crack  lengths  to  establish  a fracture  initiation  criterion  in 
terms  of  the  J-lntegral.  We  shall  show  that  there  are  some  nonlinear  visco- 
elastic materials  for  which  the  isochronal  data  In  Fig.  5 may  be  employed  In 
the  manner  described  in  [14]  for  time-independent  materials. 

It  is  sufficient  to  show  that  there  exists  a so-called  stress-working 
density  W * w(€ij>  with  the  property  that 


o 


U 


(33) 


If  W la  Independent  of  xt  the  J-lntegral  will  be  path-independent  for  two- 
dimensional  deformations  [12],  but  this  Independence  is  not  needed  for  the 
fracture  criterion.  Let  us  start  with  the  following  stress-strain  equations 
for  an  isotropic  material, 

t 

. i /■>(*  - ♦•> 


3(o*  g) 

-H— dT  * 


(34) 


(35) 


where  4 is  "reduced  time". 


F - 4(0  3 


t 


*’  3 F(T)  , 


(36) 


and  and  are  the  deviator lc  strains  and  stresses, 


eIj  5 £ij  " 3 5iJ* 


°ij  5 °lj  ‘ 3 61J* 


(37) 
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Also,  6 


1J 


Is  the  Kronecker  delta  and 


t)  £ e 


il 


9 = o 


11 


(38) 


are  the  first  strain  and  stress  Invariants,  respectively.  J(p)  Is  the  linear 
viscoelastic  creep  compliance  under  shear  stress  and  B is  the  compliance 
under  hydrostatic  loading.  Nonlinearity  Is  due  to  the  stress-dependence  of 
the  material  functions  g and  aQ.  In  view  of  the  data  on  plastics  in  [IS]  and 
the  mlcrocracklng  model  for  composites  in  [16],  It  is  assumed  that  g and  a^ 
are  functions  of  the  so-called  effective  stress 


o 

e 


1/2 


(39) 


and  the  Lebesgue  norm  of  effective  stress, 

1/q 


(40) 


These  constitutive  equations  are  special  versions  of  more  general  results 
derived  from  nonequilibrium  thermodynamics;  by  means  of  || ||  they  Incorporate 
microcracking  effects  [16].  The  equations  are  quite  realistic  In  that  they 
predict  behavior  which  is  observed  for  various  viscoelastic  materials.  Other, 
more  general  relations  could  be  employed  to  generalize  the  J-lntegral  crite- 
rion, but  we  shall  limit  our  discussion  to  the  above  ones. 

As  In  establishing  the  deformation  theory  of  plasticity,  proportional 
stressing  Is  assumed. 


°U 


■ h s 


ij  ’ 


(41) 


where  s^  Is  Independent  of  time.  Also,  let  us  assume  that  h 


ct  , where 

m • m(xi)  Is  a positive  exponent  and  c ■ clx^). 

It  la  now  possible,  with  some  amount  of  work,  to  show  that  a stress- 
working  density  exists.  Therefore,  at  any  fixed  time  prior  to  fracture  Ini- 
tiation [12] 

2 / 0 f If  « “ - s n • (42> 

o 


where  H is  specimen  thickness,  and  A ■ A(D,a,t)  Is  the  area  under  any  one 
curve  in  Fig.  6.  The  left-hand  side  of  Eq.  (42)  Is  the  J-lntegral.  Thus,  by 
observing  when  fracture  Initiates,  the  experimental  procedure  Illustrated 
above  may  be  used  to  determine  critical  J-values,  J£;  for  viscoelastic 
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Figure  6.  Force-Displacewent  Curves 
for  Various  Initial  Crack  Lengths  at 
One  Time 


materials  they  will  depend  on  the  loading  exponent  n and  on  the  time. 
It  la  important  to  observe  that  if  o{  ■ of(v), 


J 


(43) 


which  la  equal  to  2T.  Thus,  for  this  case  J is  equal  to  the  mechanical  work 
input  to  an  element  of  material  in  the  failure  zone  (per  unit  area)  at  the 
instant  of  rupture.  When  a^  ■ g ■ 1,  Eq.  (34)  reduces  to  that  for  a linear 
viecoelastlc  material,  and  it  can  be  shown  that  fracture  initiation  Eq.  (32) 
la  recovered  from  Eq.  (42)  if  Poisson's  ratio  is  one-half.  However,  when  one 
accounts  for  the  small  curvature  property  possessed  by  viscoelastic  creep 
compliancas,  as  noted  in  J5,  it  is  found  that  this  restriction  on  Poisson's 
ratio  can  be  removed  without  introducing  serious  error  in  many  cases.  Alao, 
recall  from  the  discussion  in  IS  that  Eqs.  (29)  and  (32)  often  yield  approxi- 
mately the  same  fracture  initiation  time. 


7.  CONCLUDING  REMARKS.  The  Intent  of  the  discussion  on  nonlinear  be- 
havior was  to  suggest  a procedure  for  applying  the  J-integral  to  fracture 
initiation  in  viscoelastic  materials,  and  to  provide  a theoretical  motivation 
for  this  extension.  However,  much  experimental  and  theoretical  work  is  needed 
before  one  can  have  confidence  in  the  technique.  Prediction  of  crack  speed 
in  nonlinear  viscoelastic  materials  is  a relatively  untouched  area,  except  for 
materials  in  which  viscoelastic  effects  are  limited  to  the  neighborhood  of 
the  crack  tip  (17], 

The  linear  viscoelastic  analysis  described  herein  is  limited  to  slow 
crack  growth  and  to  materials  for  which  the  neighborhood  of  the  crack  tip  is 
homogeneous.  A few  studies  exist  in  which  these  restrictions  are  not  imposed, 
but  they  are  for  very  special  cases.  For  example,  it  was  noted  that  some 
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explicit  results  exist  for  slow  crsck  propagation  between  two  different 
lncosipresslble  materials. 
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AN  INCREMENTAL  FORM  OF  THE  SINGLE- INTEGRAL  NONLINEAR  VISCOELASTIC 
THEORY  FOR  ELASTIC-PLASTIC-CREEP  FINITE  ELEMENT  ANALYSIS  * 


Duane  R.  Sanders  and  Walter  E.  Haisler 


A single-integral  nonlinear  viscoelastic  theory  is  cast  into  a generalized 
Incremental  form  for  use  in  predicting  creep  response  of  materials.  The  incre- 
mental creep  constitutive  model  is  tailored  for  use  in  an  Incremental  finite 
element  program  wherein  the  equations  of  motion  are  integrated  numerically  step- 
by-step.  Recursive  relations  for  material  memory  parameters  are  developed 
which  allow  the  retention  of  all  past  history  prior  to  the  current  time  step  and 
require  integration  over  the  current  time  step  only.  The  creep  model  parameters 
required  in  the  theory  are  obtainable  from  a standard  creep  and  recovery  test. 
The  incorporation  of  the  model  within  an  incremental  elastic-plastic  finite 
element  formulation  is  outlined. 


♦To  be  presented  at  the  3rd  ASME  National  PVP  Congress,  San  Francisco, 
California,  June  25-29,  1979. 
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INTRODUCTION 

In  a majority  of  the  finite  element  programs  in  use  today  for  analyzing  metal 
components,  the  traditional  approach  of  separating  the  total  strain  into  elastic, 
plastic,  and  creep  components  is  used  [1,2,3].  The  most  widely  accepted  means  of 
predicting  the  creep  strain  in  this  context  is  the  phenomenoloqical  creep  theory 
[4,5].  As  an  alternative  to  this  approach,  a single  integral  nonlinear  visco- 
elastic theory  derived  by  Schapery  from  thermodynamic  considerations  in  reference 
[6]  is  presented.  It  has  a form  similar  to  the  Boltzman  integral  used  in  linear 
viscoelasticity,  but  the  physical  time  is  replaced  by  a so  called  reduced  time; 
several  other  material  functions  enter  into  the  theory.  This  representation  of 
the  creep  strain  component  should  provide  within  the  same  constitutive  theory 
the  capability  of  analyzing  both  polymeric  and  metallic  materials. 

Two  constitutive  theories  similar  in  appearance  to  the  single  integral 
theory  of  Schapery' s are  given  by  Rashid  [7]  and  Valanis  [8].  Rashid  uses  a 
modified  superposition  principle  with  the  compliance,  which  is  a function  of  stress 
and  time,  transformed  such  that  it  is  a function  of  a single  parameter  called  re- 
duced time.  Valanis' s theory  [8]  is  derived  from  thermodynamic  considerations 
with  the  resulting  single  inteqral  equation  appearing  in  terms  of  the  strain 
history  instead  of  thestress  history  as  in  reference  [6].  This  theory  is  very 
similar  to  Schapery' s theory  [9]wherein  he  employed  strain  as  the  Independent 
variable  instead  of  stress.  The  functional  difference  in  references  [7]  and  [9] 
being  their  respective  definitions  of  reduced  time.  The  other  material  functions 
entering  Into  reference  [6,9]  are  not  present  in  Rashid's  or  Valanis's  theory. 

Application  of  the  single  integral  theory  has  been  made  by  Valanis  [10,11] 
for  the  viscoplastic  response  of  several  metals.  It  has  also  been  shown  to  give 
accurate  results  for  several  polymeric  composites  in  creep,  creep-recovery  and 


2 


multistep  Input  tests  on  uniaxial  specimens  [12,13]. 

A potential  advantage  of  the  single  Integral  is  that  It  easily  handles 
multistep  input  histories  including  stress  reversal  and  creep- recovery.  This  Is 

not  true  of  the  phenomenological  theory  of  creep  for  metals  based  on  strain 

i. 

hardening  which  must  resort  to  special  procedures  to  handle  the  stress  reversal 
situation  as  outlined  in  [14].  TheSe procedures  have  in  some  cases  resulted  In 
rather  poor  results  as  indicated  in  reference  [14].  Futher,  the  phenomenological 
theory  has  no  predictive  capability  for  creep-recovery  input  stress  history.  Thus, 
the  single  integral  constitutive  theory  may  offer  aaditonal  predictive  capabilities 
over  the  phenomenological  theory  of  creep. 

To  put  the  single  integral  theory  in  a form  suitable  for  use  with  a incre- 
menatal  finite  element  program  wherein  the  equations  of  motion  are  integrated 
..  numerically  step-by-step,  the  time  derivative  of  the  single  integral  equation  is 

taken.  In  this  form,  the  creep  strain  increment  can  be  determined  for  the 
current  time  interval  based  on  the  creep  strain  rate  at  the  beginning  of  the  step. 
Computationally,  the  relationships  derived  for  the  rate  form  of  the  single  integral 

i 

equation  are  particularly  appealing  since  recursive  relationships  are  developed 

1 

which  allow  for  the  retention  of  all  past  history  prior  to  the  current  time  step 
and  require  integration  over  the  current  time  step  only. 

The  incorporation  of  the  model  within  an  incremental  elastic-plastic  finite 
element  formulation  is  outlined.  The  time-independent  Incremental  plasticity 
relations  with  a Von  Mises  yielo  criterion  and  a kinematic  hardening  rule  are 
employed  to  describe  the  plastic  component  of  strain. 


A graphical  procedure  for  determining  the  material  functions  from 
either  creep  data  or  creep  equations  is  presented.  It  is  shown  that 
this  procedure  is  only  applicable  to  materials  whose  elastic  strain  on 
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unloading  Is  equal  to  that  on  loading. 

Numerical  comparison  of  the  single  integral  theory  to  that  of  phenmenological 
theory  of  creep  are  given  for  creep,  creep  recovery,  multi  step  and  reverse  loading 

histories  for  a type  304  stainless  steel.  The  numerical  results  were  obtained 
using  the  computer  program  NONVIS  which  implements  the  integration  procedure 

developed  below  for  the  single  integral  creep  equation. 

It  is  suggested  that  the  rate  from  of  the  sinqle-integral  equations 

appearing  in  [7,8,9]  could  be  integrated  by  employing  a procedure  similar  to 
that  given  in  this  paper  for  the  single  integral  constitutive  theory  of  references 
[6],  provided  the  kernel  function  is  a combination  of  constants,  linear  , and 
exponential  terms  as  will  be  outlined  below. 


RATE  FORM  OF  THE  SINGLE  INTEGRAL  EQUATION 

The  single  integral  nonlinear  viscoelastic  constitutive  equation  derived 
in  reference  [6]  for  an  initially  isotropic  material  is  written  for  isothermal 
conditons,  incompressible  creep  strains,  and  infinitesimal  strains  as. 
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and  A^,  and  xr  are  constants.  The  functions  ag,aDand  oC(j  are  material 

properties  and  are  in  general  function  of  the  three  stress  invariants.  eab  are 
the  engineering  strain  components  and  are  the  components  of  the  Cauchy  stress 
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tensor.  The  function  AJU)  is  the  transient  component  of  the  so  called  linear 
viscoelastic  shear  creep  compliance.  The  first  term  in  equation  (1)  represents 
the  instantaneous  elastic  response  of  the  material. 

To  cast  equation  (2)  into  a form  suitable  for  use  Ip  a incremental  finite 
element  analysis  the  time  derivative  of  equation  (1)  is  taken,  and  equation  (2a) 
is  used  to  obtain. 


to  ODtain,  a t 
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Equation  (3)  is  rewritten  as, 


f b^cJi  \ 

= T(o)  ^ Xj  , 
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where. 


at 
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and 
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A J'lt-VJ  * — 
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A time  marching  procedure  for  evaluating  X ^ and  YC{j  is  formulated  below. 
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Evaluation  of  Xcd  and  Ycd 


The  evaluation  of  the  two  integrals  Xcd  and  Ycd  given  by  equations  (5) 
and  (6)  are  put  in  a form  suitable  for  use  in  an  incremental  finite  eianent 
analysis.  Once  the  value  of  the  integrals  are  known,  the  creep  strain  rate 
is  evaluated  from  the  sum  of  the  second  and  third  terms  of  equation  (4)  as: 

/ ^ \ A 

arson:/  v . . do  y 

• Tt  X‘J  + te;  • (8) 


Two  important  features  in  the  evaluation  Xcd  and  Ycd  are  the  determination  of 

two  recursive  formula's  whose  values  are  known  from  the  previous  time  step  and 

the  evaluation  of  the  two  integrals  which  require  integration  only  over  the  current 

time  step.  The  formulation  presented  is  similar  to  the  integration  technique 

employed  by  Zak  and  Hermann  in  reference  [15,16]  for  a thermorheological  simple 

material  due  to  a variable  temperature  field. 

Subdivide  the  time  from  0 to  t into  N subintervals  such  that 

t0  = = t0  + 4t],  . . .,  tN  =tN-1  + At  = t.  For  ocd(x,tN)  andld(x,tN) 

write  a . and  and  assume  that  (0,tM)  = 0 and 
cdN  cdN  caN  H 

ecd(0,t  ) = 0.  Then  X . and  Y . can  be  written  as 
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Now,  rewrite  equation  (13)  using  equation  (16),  as, 
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Where, 
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Recursive  relationships  for  equations  (23)  and  (24)  are  developed  In 
Appendix  I and  are  given  below  as: 
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where  , <'2'c^w.l  , ^iw.l  a*A  TZ  are  known  from  the  previous 

Tt(r> 

time  step.  Thus,  equation  (19)  can  be  evaluated  by  integrating  and 

-r?  (r,) 

equations  (17)  and  (18)  respectively,  over  the  current  time  step 

only.  All  history  effects  are  contained  in  L-j  .equation  (20).  In  Appendix  I, 

N 

Equations  (17)  and  (18)  are  set  up  for  numerical  integration  using  Gaussian 
Quadrature.  The  resulting  equation  for  equation  (17)  is: 
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where. 
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and 


The  resulting  equation  for  equation  (18)  is: 


72 


= % - J2 


(1.20) 


where. 


J ^ + t"  +t^-')  A*  (1.19) 
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Evaluation  of  the  reduced  time  integral  equation  (2a)  is  also  presented  in 
Appendix  I. 
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Evaluation  of  Yccj 

N 

Following  a procedure  similar  to  that  used  in  evaluating  XC{j  above,  the 
intergration  of  YC(^  is  performed.  YC(^  is  given  by  equation  (10)  as. 


where. 


A,r>  -VXr  „ (r.j 

-ZB 

r r' 


(10) 


(25) 


Expanding  equation  (10)  as  the  sum  of  two  integrals,  using  equation  (25),  and 
approximating  the  time  rate  of  change  of  GC(j  given  by  equation  (12),  equation  (10) 
is  rewritten  as, 
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Defining  the  following  quanities, 
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Note  that  JW  ' is  the  same  quantity  that  appears  in  the  evaluation  of 
Xc . . Using  equation  (27)  and  (28),  equation  (26)  is  rewritten  as. 


U - ^ 7"  + 


A recursive  relationship  for  equation  (31)  is  developed  in  Appendix  I and 
is  given  below  as: 


F.  * e 


rlr'  Titr) 

where  J n-I  are  '<nown  ^rom  the  Previous  time  step.  Thus  Yc^  . 

equation  (31),  can  be  evaluated  by  intergrating  Jl^,  equation  (17)  over 

the  present  time  step  only.  All  history  effects  are  represented  by  K . , 

N 

equation  (28).  In  Appendix  I,  equation  (17)  is  set  up  for  numerical 
integration  using  Gaussian  Quadrature  (the  results  were  given  above  by 
equations  (1.15)  and  (1.16). 
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INCREMENTAL  EQUATIONS  OF  EQUILIBRIUM 


In  this  section,  an  incremental  finite  element  formulation  for  the 
isothermal,  elastic-plastic-creep-large  strain  problem  is  developed.  The 
formulation  is  suitable  for  use  with  any  creep  constitutive  theory  which 
will  predict  the  incremental  creep  strain  component  based  on  the  creep  strain 
rate  the  beginning  of  the  step.  Although  the  creep  and  plasticity  relations 
are  for  small  strains,  the  formulation  is  developed  for  large  strain.  One 
may  begin  with  the  equations  of  equilibrium  written  in  terms  of  the  second 
Piola-Kirchoff  stress: 


fa + 2*.*)]  f f*  = ° 


where  aj  and  u^  are  Lagrangian  coordinates  and  displacements,  respectively, 
S^is  the  2nd  Piola-Kirchhoff  stress  tensor,  pQ  is  undeformed  density  and 
is  the  body  force  per  unit  undeformed  volume  per  unit  mass.  Applying  the 


virtual  work  principle  at  time  t + At  yields 


t+  At 


where 


t *4  At  t+At  < 

1, 


+4tt  I t+*t  t4it  f t+At  t+** 

JR  rjsT;  Jkt  AS  + ^ f.  ^ M ( 

and  6E*jAt  Is  the  variation  of  the  Green-Lagrange  strains  at  time  t+At  and 
T*+At'are  surface  tractions  at  time  t+At  applied  to  the  deformed  surface  S. 
Equation  (2)  may  be  put  into  Incremental  form  by  writing 
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where  and  are  stresses  and  strains  at  time  t and  aS^  and  AE^j  are 
increments  of  stress  and  strain,  respectively.  The  strain  increment  may  be 
decomposed  into  components  which  are  linear  and  nonlinear  in  the  displacement 


increments 


AEv  - * A£< 


Where  ^ £AW,  * ^ v 

wc  | ( ^ 

AE-  r 2.V  5 A.-  / 

Substituting  equations  (35)  and  (36)  into  (33)  yields 

f * < , f ^ . L tlL  tf&  t 

iVV'  + 1 % ,* 

vp  % *4 

The  stress  increment  may  be  decomposed  into  two  components,  one  which  Is 
dependent  upon  total  strain  and  one  which  is  independent  of  strain  (i.e., 
creep , thermal , etc . ) : 

- V + ^ <3! 

where  D.j^ls  the  usual  effective  tangent  modulus  and  aP^  Is  a stress 
Increment  due  to  strain  independent  phenomena  (as  Is  usually  assumed  in 
creep).  Substituting  equation  (39)  into  equation  (38),  making  use  of 
equation  (36),  and  neglecting  terms  which  would  be  nonlinear  In  displacement 
Increments,  yields  the  following: 


i. 


13 


. - 1 (U)  * aPo)  'K)"  + st 


t+^t 


The  term  aP^j  may  be  interpreted  as  the  change  in  stress  required  to 
account  for  the  creep  and  thermal  strains.  Equation  (40)  takes  on  the 
following  form  when  put  into  matrix  form 

[*] {**'“]  - ( [k‘J  * 1 *SJ)W1  (41) 

where  [M]  is  the  mass  matrix,  [K*]  and  [K^]  are  "linear"  and  "nonlinear" 
stiffness  matrices,  {Rt+At}  is  a vector  of  forces  due  to  externally  applied 
loads,  {F*}  is  a vector  of  forces  due  to  internal  stress,  and  { Aq>  is  the 
Increment  of  the  nodal  displacements  *.  Complete  details  of  the  derivation 
of  the  quantities  in  equation  (41),  without  creep,  may  be  found  in  reference  £17]. 

We  now  present  a summary  of  the  determination  of  and  AP^  for 

kinematic  hardening  for  the  isothermal  case,  we  assume  a yield  function  can  be 


expressed  by 


F-  i - K 


1 - o 


For  kinematic  hardening,  write 


•f  ■ i ( V % ] ( s.v ' < > 


and  Sjj  and  are  devlataric  components  of  the  2nd  Piola-Kirchoff  stress  and 
yield  surface  center,  respectively  and  K is  a constant. 

The  associated  flow  rule  Is  then  given  by 

jF  <*♦> 
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where  X is  a scalar  to  be  determined.  Using  the  Ziegler  modification  of  the 
Pr  ger  kinematic  workhardening  rule,  it  is  assumed  that  the  projection  of  the 
increment  of  stress  onto  the  normal  to  the  yield  surface  in  a scalar  multiple 
of  the  dot  product  of  the  plastic  strain  and  yield  surface  normal 

ar  ,rf  iF 

<)%  = d dt‘l  iSij  (45) 


where  C is  a scalar  (hardening  modulus  ) to  be  determined  from  a uni-axial  stress- 
strain  curve.  It  can  be  shown  that  equation  (45)  is  equivalent  to  the 
consistency  condition  for  kinematic  hardening  and  isothermal  conditions.  For 
small  strains,  the  decomposition  of  strain  is  assumed: 


ds-  = 


where  is  the  elastic  constitutive  tensor  and  dE  £ are  the  Incremental 

Ijmn  mn 

creep  strain  components.  The  incremental  creep  strain  components  are  deter- 
mined from  equation  (8)  and  a time  increment  over  which  it  is  assumed  the 
creep  strain  components  are  constant.  Substituting  equation  (44)  into  equation 
(45)  and  (46)  and  solving  the  resulting  equation  for  x, 
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Substituting  equation  (47)  into  equation  (44)  and  the  results  into  equation 
(46)  and  comparing  the  results  with  equation  (39)  one  can  show  that  the 
Instantaneous  modulus  tensor  is  given  by 
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and  the  change  of  stress  due  to  creep  strains  is  given  by 


^ ■ 


Equation  (46)  can  now  be  written  as. 
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DETERMINATION  OF  MATERIAL  FUNCTIONS  FROM  EXPERIMENTAL  DATA 

In  this  section,  equation  (1)  is  specialized  to  the  case  of  uniaxial 
loading;  and  to  materials  whose  elastic  response  on  loading  and  unloading  in 
a creep  and  recovery  test  are  the  same.  The  required  material  functions  can 
then  be  evaluated  by  a graphical  procedure  from  the  results  of  creep  tests. 

For  a general  description  of  the  methods  for  obtaining  the  material  functions 
from  experimental  data,  see  Lou  and  Schapery  [13]. 

Consider  a uniaxial  bar  which  is  in  a uniform  uniaxial  state  of  stress  at 
a uniform  temperature  considering  only  one  stress  o,  and  one  strain  e . we  let 
e=  e 1 1 and  a * a^.  When  equation  (l)  is  specialized  to  this  uniaxial  loading, 
the  creep  term  is  given  by 


and  the  three  Independent  material  functions  are  given  by 
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It  is  shown  In  reference  [12]  that  g^a)  = 1,  implies  that  the  amount  of 
elastic  straining  on  unloading  is  the  same  as  on  loading  and  that  ®aW  * 

A 

Using  g,(rf)  = 1 and  ? <T4b  , 


equation  (51)  is  now  given  by, 
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where 


Ur)  s 


<T 

a,. 


and 


(57) 

(58) 


G.6(<r) 

When  equation  (56)  is  subjected  to  a creep  test  the  resulting  equation 
can  be  used  to  evaluate  the  material  functions  g2  and  ao.  Consider  equation 

(56)  subject  to  the  stress  history: 

<r(t)  - H (t)  <r  . (59) 

Substituting  equation  (59)  into  equation  (56)  gives 


£c  * ^ (X-l  r ; 

and  the  creep  compliance  is  given  by 


(60) 
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Taking  the  log  of  both  sides  of  equation  (61)  gives. 


ad(*)  - 


The  simple  form  of  equation  (62)  is  such  that  a simple  graphical  procedure  can 
be  employed  to  evaluate  the  material  functions.  This  is  easily  demonstrated 
by  plotting  AD  and  t on  double  - logarithmic  paper,  then  curves  at  different 
stress  levels  can  be  superimposed  by  translating  them  along  the  AD  and  t axes. 
If  the  curves  are  shifted  to  a reference  stress,  , the  amount  of 

horizontal  shift  (t)  and  vertical  shift  (aD)  equals  log  g?  = * —7—.  and 

logar=  log  ^ , respectively.  These  shifts  are  graphically  demonstrated 

in  figure  ( 1 ). 


The  direction  at  the  shifts 
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are  as  follows: 
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SOME  NUMERICAL  RESULTS 


(63) 


This  section  presents  some  numerical  results  obtained  for  a type  304 
stainless  steel  with  the  computer  program  NONVIS  which  implements  the  pro- 
cedures developed  above  for  integrating  equation  (1)  for  the  following  stress 
histories:  creep  (figure  2),  multi-step  (figure  3),  creep-recovery  (figure  4), 
and  stress  reversal  (figure  5).  The  results  obtained  are  compared  with  those 
of  the  phenomenological  creep  theory  based  on  strain  hardening  using  total 
creep  strain  [18]. 

An  effective  creep  equation  of  the  form: 

--  A (/-  ) + Kt  (64) 


was  used  to  determine  the  material  functions  in  lieu  of  the  actual  experimental 
data.  Two  different  heats  of  type  304  stainless  steel  are  considered.  The 
creep  strain  parameters  for  heat  no.  9T2796  at  1100°F  are: 
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with  <r  (effective  stress)  in  ksz,  and  t in  hours.  The  creep  strain  parameters 


for  hear  No.  8043813  at  1200UF  are: 


and 
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The  vertical  shift  is  defined  by  the  function  A(r)  and  the  horizontal  shift 
by  the  function  r(f).  This  can  be  seen  by  equating  the  transient  term  given  in  - 
equation  (60)  (for  r = 1,  r'=  0)  and  (64)  from  which  it  is  determined  that 


4*  Sx  r 

(67a) 
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(67c) 

and  ?r  s 1-0 

(67d) 

Rearranging  equation  (67a)  gives 

' V - , 

(68) 

f 

and  from  equation  (2a)  and  (67b) 
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it 
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Therefore,  the  actual  shifting  of  the  creep  compliance  wasn't  necessary  to 

determine  the  material  functions  g 2(<r)  and  £f).  it  is  easy  to  check  the 

function  g and  given  above  by  using  the  graphical  procedure. 

2 

If  the  steady  state  creep  term  in  equation  (65),  term  linear  in't  , is 
subjected  to  the  creep  stress  history,  equation  (59),  it  gives, 

e<°  • <70> 

Using  equations  (68)  and  (69),  equation  (70)  becomes 

f = Art  (71) 

cs 

Equation  (64)  above  gives  for  the  steady  state  creep  term 


f,  ' <t.  (72) 

However,  from  equations  (71),  (72),  and  (65c)  it  is  seen  that  equations 
(71)  and  (72)  are  not  equal.  Thus  it  is  concluded  that  the  reduced  time 
cannot  be  related  through  the  reduced  time  use  for  transient  creep. 

In  the  application  to  the  type  304  stainless  steels,  the  steady  state 
term  will  be  given  by  the  second  term  of  equation  (64),  £jsK.t. 

Figure  2 presents  creep  strain  vs.  time  for  a uniaxial  stainless  steel 
(heat  No.  9T2796)  specimen  for  a creep  test.  The  results  of  both  the 
nonlinear  viscoelastic  and  the  phenomenological  creep  theory  are  nearly 
identical  as  is  expected.  Only  the  nonlinear  viscoelastic  theory  is  presented. 

Figure  3 presents  the  creep  strain  response  of  a type  304  stainless 
steel  (heat  No.  8043813)  subjected  to  a creep  and  recovery  tests  for  a 


Creep  Strain  (in/ in)  Creep  Strain  (in/in) 
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Figure  2.  Creep  of  A Uniaxial  Specimen 


Figure  3.  Creep-Recovery 
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uniaxial  specimen.  The  results  Indicate  very  little  recovery;  of  course,  the 
phenomenological  theory  predict  no  recovery.  However,  It  should  be  noted  that 
the  nonlinear  viscoelastic  theory  would  predict  complete  recovery  If  the 
creep-recovery  curve  were  extended  to  larger  values  of  time. 

Figure  4 presents  the  creep  strain  response  of  a type  304  stainless  steel 
(heat  No.  8043813)  uniaxial  speciman  subjected  to  montonlcally  Increasing  load. 

The  results  are  In  reasonable  agreement  with  the  experimental  results  and 
are  very  close  to  the  results  base  on  strain  hardening. 

Figure  5 presents  the  creep  strain  response  of  a uniaxial  specimen  subjected 
to  a reversed  cyclic  loading  history.  On  the  first  reversal  the  nonlinear 
viscoelastic  theory  over  predicts  the  creep  response  as  given  by  the  experimental 
data.  This  is  due  in  part  to  the  lack  of  significant  recovery  for  Type  304  stain- 
less steel.  To  study  the  behavior  consider  the  stress  history  up  to  the  time 
of  the  second  stress  reversal , 

r(i)--rH(t)  - **(*"*•)  , o<t<tv  (73) 

Substituting  equation  (73)  into  equation  (56)  gives 

= D(*)  gx«r)  - c^v(<r,)  D(*-V.)  . (74) 

Thus  it  is  seen  that  the  second  term  gives  2g2( o i ) times  the  compliance.  Thus 
twice  the  value  of  g2(o])  is  multiplying  the  compliance  than  was  on  initial 
loading.  This  seems  to  be  okay  for  materials  with  an  active  recovery,  l.e., 
polymers.  However,  for  metals  the  recovery  is  only  slight  and  the  second  term 

over  predicts  its  contribution  to  the  recovery.  Futher  the  first  term  is  related 
to  the  concept  of  a fading  memory,  which  is  questionable  for  the  stress  reversal 
situation.  Note  the  curve  plotted  for  the  first  stress  reversal  with  no 
recovery  or  memory  of  past  history.  This  curve  was  obtained  by  initializing 
equation  (20)  to  zero,  which  represents  the  memory  of  the  material  to  prior 


- ■ .. 


Creep  Strain  (in/ in) 


24 


0 200  400  600  800 

Time  (hr$) 

Figure  5.  Creep  Strain  Response  of  Uniaxial  Specimen 
Subjected  to  Cyclic  Load. 


histories,  and  setting  G„  =0  when  a stress  reversal  occurs. 

N-2 

same  response  as  an  initial  loading  with  the  strain  accumulated 


This  gives  the 
at  the  time 


of  the  reversal  serving  as  a new  origin. 


CLOSURE 


The  purpose  of  this  paper  was  first  to  present  an  efficient  method  of 
Integrating  Schapery's  single  integral  equation  for  use  In  a incremental  finite 
element  program  wherein  the  equations  of  motion  are  integrated  step-by-step. 

And  secondly  to  present  some  preliminary  results  for  metals  using  this  theory. 

It  is  believed  that  the  method  presented  for  integrating  the  single  Integral 
equation  is  efficient  based  on  the  fact  that  the  integrations  required  are  for  the 
current  time  step  only  and  does  not  require  integration  over  the  entire  past 
history  to  determine  the  current  creep  rate.  From  the  experience  gained  with 

NONVIS,  the  programming  aspect  of  implementing  the  theory  is  straight 
forward  and  can  easily  be  incorporated  into  existing  finite  element 
programs . 

The  results  of  this  paper  are  not  presented  to  draw  final  conclusions 
regarding  whether  Schapery's  nonlinear  viscoelastic  theory  are  the  phenomenological 
theory  based  on  strain  hardening  is  superior.  Obviously,  the  application  of  the 
phenomenological  theory  for  metals  has  developed  over  many  years  and  the  application 
of  nonlinear  viscoelastic  theory  is  just  beginning  its  development  for  metals. 
However,  from  the  results  presented  here,  it  is  seen  that  in  all  except  the  case 
of  stress  reversal  the  nonlinear  viscoelastic  theory  gives  results  as  acceptable 
as  those  based  on  strain  hardening  theory  and  in  fair  agreement  with  experimental 
data.  For  the  stress  reversal  case  two  problems  arise.  The  nonlinear  viscoelastic 
theory  predicts  a complete  recovery  of  the  (as  in  linear  viscoelasticity)  strain 
accumulated  at  the  first  and  subsequent  reversals.  Secondly,  the  theory  has  a 
fading  memory  of  past  stress  history.  However,  as  was  noted  these  effects  can  be 
neglected  for  the  stress  reversal  case  by  initializing  to  zero  the  parameters 
that  contain  the  memory  of  past  histories  when  a stress  reversal  occurs  and 

setting  the  value  of  G,,  =0. 
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Although  the  results  obtained  have  been  for  the  response  of  a stainless 
steel  most  of  the  results  reported  to  date  have  been  for  a polymer  or  polymer 
composite.  As  noted  In  the  Introduction,  reference  [12,13]  show  the  theory  In 
good  agreement  for  a polymer  composite  for  the  following  stress  histories:  creep, 
creep-recovery,  and  multistep  inputs.  Therefore,  it  is  believed  that  the  use 

. . 

of  the  integration  technique  for  the  single  integral  equation  in  conjunction  with 
a incremental  finite  element  analysis  has  immediate  application. 

Current  work  with  the  single  integral  equation  of  Schapery’s  involves:  (1) 
providing  a better  predictive  capability  for  the  stress  reversal  and  other 
stress  histories  for  metals,  (2)  determining  the  creep  compliance  and  material 
function  for  several  composite  materials  and  (3)  incorporating  the  above  Integra- 

I | T 

tion  procedure  into  a finite  element  program  for  testing  on  two  dimensional 

* » 

problems. 
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APPENDIX  I 

«■>  -(r') 

Recursive  Relationships  for  Clcd*  and  CZcJ* 

After  expanding  equation  (23)  as  the  sum  of  two  integral,  it  can  be 


rewritting  as 

Cl 


<r>  _ -(vV,)/rr  ( (tw*1  -<wvyr,  AhLjt' 

- e U e 
+ *”•*'>/*  ^ 


(1.1) 


To  evaluate  equation  (1.1)  note  that  if  in  equation  (23)  N is  replaced  by  N-l, 


equation  (23)  can  be  written  as. 


which  is  the  value  of  Cl  from  the  previous  time  step. 


(1.2) 


Define  as. 


(1.3) 


and  approximate  the  time  rate  of  change  of  Gc(Jas 

at  * 4 V, 

on  the  time  interval  t N _ £<  * < ^-T  ^sing  equation  (1.2),  (1.3)  and 
(1.4)  in  equation  (1.1),  gives. 


(1.4) 


ci,  = e 

ctH~\ 


tr>  + [ e 

A Xff  x 


In  the  integral  given  by  equation  (17),  replace  N by  N-l  to  get, 

n*  ■ t r . 


* (1.5) 


(1.6) 


N-l  yi 


I 
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Which  is  Jl^rVrom  the  previous  step.  Thus,  equation  (1.5)  can  be  written  as, 

(r ) 

■ e i ci 

VN 

,ir) 


(r) 


- i^'r'  { cC"  + KVl  + W)  «r„„]  (I. 

Note  that  the  values  of  Cl(j 


7) 


AM 


^ , 6cj  and  ^ are  known  values 

A l-l  V'J  N-2. 


from  previous  steps. 

Similarly,  equation  (24)  can  be  manipulated  and  rewritten  as, 

C2.  = ■+  CZ.j  f-  ( 6.  i *6  . ) T<2r<  (1.8) 

'"'Sd*  *V-»  c<Vi  1 cVi  r*-i  ' ' 

(f‘)  xi 

Where  CZ( ^ , and  are  values  of  equation  (24)  and  (18)  from  the 


previous  step. 


Recursive  Relationship  for  F Cci 


( r ) 


After  expanding  equation  (30)  as  the  sum  of  two  integrals.  It  can  be  written  as, 
(r>  _ rtv-t  -(*„  -V*_,  r V,  ♦ 46^  J , 

'*4*  “ J0  dv 


+ r4"'1  ■(*»-%-,  -fft  ^ 


■J 


(1.9) 


To  evaluate  equation  (1.9)  note  that  if  in  equation  (30)  N is  replaced  by  N-l, 
equattion  (30)  can  be  rewritten  as, 

tr(r)  (e*1  e ^*’’  *'^tr  j6cJ  dt 


a 


4t' 


AM 


(I. 10) 


which  Is  the  value  of  from  the  previous  step.  Using  equations  (1.9). 

(I. 10)  and  (1.4),  equation  (30)  can  be  written  as,  . 

. . . . _ - / -I 


ir)  -*+*/*,{  ptr'  + ~ (1.11) 

^ r e l J"-'  + At»l  Jtf.i 

The  Integral  term  in  equation  (1.14)  can  be  written  as  Jl^j  using  equation 


(1.6)  thus 


< * <]  • 


1 


(1.12) 
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(r)  (r 1 ) 

Evaluation  of  the  quanities  J1  and  J2 
" N N 


,<r) 


Expanding  J1  N , equation  (17)  as 

U * •" 


Jlr.  - «, 


or 


r At#/ 


f*'  e 

*-  ... 

" Tr  At 


|C 


ty  r -(V*'Vrr  - 


♦n« 

H-\ 


- e 


» 


(1.13) 


Let  the  change  in  reduced  time  during  the  interval  At^  be  related  to  At ^ through 
the  parameter,  a 


% - Vf  = ^ or 

A - 

AV/r 

(1.14) 

Define  the  quantity  as 

. ffc*r  -(V 

TI*'rU  \ c - e 

-t')] 

at' 

(1.15) 

Thus,  ,/>  JL  £ f + 

Jl  = At*.  O-  1 ” J i 

< 

(1.16) 

M * L t/*-* 

Tr  -AVv/r  , Ti* 

* 7f„  (<-  e > ¥ Jir- 

Make  the  following  change  of  variable  in  J1  N 


t - 


t-U  i t»i 


A*v 
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at  lower  limit 


at  upper  limit 


t * *f 


2 


or  Jt  - ~ r 4 * 


4ent< 


2t  - - **•' 

4 tv 

« r 6**+**  * *«i' 


Jif  - i j .e-v*. 


r 4 tv  - tv  - 1 


‘v  " ‘v-;  i , 

* ' J It 


(1.17) 


Equation  (1.17)  is  set  up  so  Gaussian  Quadrature  can  be  used  for  numerical 
evaluation. 

Expanding  equation  (18)  as 


(1.18) 


■ iy jt  ■ 
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r1* 

n//.  iA  + At1 


(1.23) 


l*-l 


at  a given  point  in  space  approximate  the  stress  as  a linear  function  of  time 


with  the  Interval  At„ 


o s ' 

N-l  ^ y I 

A ''  ! 

u I 


A1N  I 

. » 


t t, 


Thus, 


<T(t)=  A,  -4-  A,  t 


i t 5 


where 


. <Vi  *»>  “ *nji  A,  ; 


At,» 


change  the  variable  of  integration  to  o,  and  obtain 


f 8 1- 


, & 

Mr'>  Al- 


(1.29) 


